panel.fit: Panel functions for predicted values and SE bands

View source: R/gpanel.R

panel.fitR Documentation

Panel functions for predicted values and SE bands

Description

Panel functions for predicted values and SE bands using 'layer' and 'glayer' in the package latticeExtra

This is an experiment in writing a function that can be called via layer or glayer without further complications e.g. xyplot(......,labels = rownames(data)) + layer(panel.labels(...)) or xyplot(....., labels = rownames(data), subscripts = T) + glayer(panel.labels(...)). For selected labels see the examples with trellis.focus and panel.identify

Fills in values in gaps between observed predictor values to help produce a smooth graph of predicted values versus predictor values.

Usage

panel.fit(
  x,
  y,
  fit,
  lower,
  upper,
  subscripts,
  ...,
  type,
  group.number,
  alpha,
  col,
  col.line,
  col.symbol,
  border = F,
  font,
  fontface
)

gpanel.fit(
  x,
  y,
  fit,
  lower,
  upper,
  subscripts,
  ...,
  type,
  group.number,
  alpha,
  col,
  col.line,
  col.symbol,
  border = F,
  font,
  fontface
)

panel.band(
  x,
  y,
  fit,
  lower,
  upper,
  subscripts,
  ...,
  type,
  group.number,
  alpha,
  col,
  col.line,
  col.symbol,
  border = F,
  font,
  fontface
)

panel.errorbars(
  x,
  y,
  fit,
  lower,
  upper,
  length = 0.2,
  angle = 90,
  subscripts,
  ...,
  type,
  group.number,
  alpha,
  col,
  col.line,
  col.symbol,
  border = F,
  font,
  fontface
)

panel.labels(x, y, labels, subscripts, ...)

fillin(data, form, n = 200, xpd = 1)

Arguments

x, y

position of labels, usually supplied through panel call

y

normally passed by 'layer' or 'glayer'. Used only in 'panel.labels'.

fit

fitted values of a model, generally passed through 'layer' from a call to 'xyplot': e.g. xyplot( y ~ x, data, groups = g, fit = data$yhat, lower = with(data, yhat - 2*se), upper = with(data, yhat + 2*se), subscripts = T)

lower, upper

lower and upper limits of error bands, passed from main plotting function

subscripts

subscripts, passed from main plotting function

...

NOTE: may specify anything you don't want passed through ...

group.number

, passed from main plotting function

alpha

transparency, passed from main plotting function

col

passed from main plotting function

col.symbol

is used to control color when using 'groups'

border

default = FALSE for panel.band.

font

passed from main plotting function

fontface

passed from main plotting function

length

(default .2) of crossbar in panel.errorbars

angle

(default 90) of crossbar in panel.errorbars

labels

to display

data

data frame with values of x that need filling in

form

formula identifying variable x to fill in and grouping variables, g1, g2, etc. using syntax: ~ x + g1 + g2 (the variable to fill in comes first)

n

number of additional points over range of predictor (default 200)

xpd

expansion beyond range of predictor (default 1.0, i.e. no expansion)

X

argument, passed by 'layer' or 'glayer'

dots

passed from main plotting function

Details

With 'layer' and 'glayer' in 'latticeExtra', these functions can be used to easily generate fitted values and confidence or prediction bands that have a reasonable appearance whether a plot uses 'groups' or not.

Value

The 'panel.bands', 'panel.fit', 'panel.errorbars', and 'panel.labels' functions are invoked for their graphical effect.

Functions

  • gpanel.fit(): to be used with 'layer' – but, actually, identical to 'glayer'

  • panel.band(): identical to panel.fit but kept for backward compatibility

  • panel.errorbars(): similar to panel.fit but draws error bars instead of bands

  • panel.labels(): add labels

  • fillin(): fill in values to produce smoother fitted curve

Author(s)

Georges Monette georges@yorku.ca

Examples

## Not run: 
  library(spida2)
  library(latticeExtra)
  library(car)
  
  ### Show data, predicted values and confidence bands 
  ### using a 'predict' method that produces 'se'   
  ### -- showing predicted values only at points in the data

  fit <- lm(prestige ~ (income + I(income^2)) * type, Prestige,
      na.action = na.exclude)
  pred <- cbind(Prestige, predict(fit, newdata = Prestige, se = TRUE))
  head(pred)
  (p <- xyplot( prestige ~ income , pred,  groups = type,
                subscripts = T,
                fit = pred$fit,
                lower = with(pred, fit - 2*se.fit),
                upper = with(pred, fit + 2*se.fit)))
  p + glayer(panel.fit(...))

  ###
  ### Use 'fillin' to add points in sparse regions of the predictor
  ### to produce smoother bands
  ###

  fit <- lm( income ~
        (education+I(education^2)+I(education^3) +I(education^4))* type,
        Prestige, na.action = na.exclude)    # overfitting!
  # adding extra values of predictor to get smooth line
  Prestige$occupation <- rownames(Prestige)
  z <- fillin(Prestige, education ~ type, xpd = 1.1) # fill in 'education' within 'type'

  dim(z)
  dim(Prestige)
  z <- cbind(z, predict(fit, newdata = z, se = TRUE))
  head(z)
  gd(3, cex = 2, lwd = 2, alpha = .7)
  (p <-  xyplot( income ~ education, z, groups = type,
                 subscripts = T,
                 fit = z$fit,
                 lower = z$fit - z$se,
                 upper = z$fit + z$se,
                 auto.key = list(space='right', lines= T)))
  p + glayer( panel.fit(...))
  p + glayer( panel.fit(...,  alpha = .1))
  # Using spida2:gd() to get a ggplot-like appearance
  gd(3,lty=1,lwd=2)
  p + glayer( panel.fit(...,alpha = .1))

  ###
  ###  Prediction with 'wald' for fits with predict
  ###  methods that don't produce se's, e.g. lme
  ###
  
  library(nlme)
  fit <- lme(mathach ~ (ses+I(ses^2)) * Sex * Sector, hs, random = ~ 1|school)
  summary(fit)
  pred <- expand.grid( ses = seq(-2,2,.1), Sex = levels(hs$Sex), Sector = levels(hs$Sector))
  head(pred)
  # merge pred with original data
  hs$from <- 'data'
  pred$from <- 'pred'
  dm <- merge(hs, pred, all = T)
  w <- wald(fit, getX(fit, data = dm)) # attaches data to wald.object so it can included in data frame
  # or:
  w <- wald(fit, pred = dm)  
  w <- as.data.frame(w)
  head(w)
  library(latticeExtra)
  gd(pch = 1, alpha = .2)
  (p <- xyplot(mathach ~ ses | Sector, w, groups = Sex,
       auto.key = T, type = 'p',
       fit = w$coef,
       upper = with(w,coef+2*se),
       lower = with(w,coef-2*se),
       subscript = T))
  p + glayer( gpanel.fit(...,alpha=1))

  wald(fit, 'Sex')  # sig. overall effect of Sex
  wald( fit, ':Sex') # but no evidence of interaction with ses
  wald( fit, '\\^2') # nor of curvature

  # but we continue for the sake of illustration

  L <- Lform( fit, list( 0, 1, 2*ses, 0, Sex == 'Male', (Sex == 'Male')*2*ses), hs)
  L
  (ww <- wald ( fit, L ))
  wald.dd <- as.data.frame(ww, se = 2)
  head( wald.dd )

  xyplot(coef ~ ses | Sex, wald.dd, type = 'n',
      ylim = c(-5, 12),
      ylab = 'increase in mathach per unit increase in ses',
      fit = wald.dd$coef,
      upper = wald.dd$U2,
      lower = wald.dd$L2,
      subscripts = T) +
  layer(panel.fit(...))


  ###
  ###  Using panel.fit with no groups
  ###

  (p <-  xyplot( income ~ education| type, z,
                 subscripts = T,
                 fit = z$fit,
                 lower = z$fit - z$se,
                 upper = z$fit + z$se,
                 auto.key = list(space='right', lines= T)))
  p + layer( panel.fit(...))

  gd_(basecol = 'tomato4')  # Use 'gd_' to set parameters without groups
  p + layer( panel.fit(...))
  p + layer( panel.fit(...,  col = 'grey10'))

  ###
  ### With panels and groups
  ###

  z <- Prestige
  z$gender <- with(z, cut( women, c(-1,15,50,101),labels = c("Male","Mixed","Female")))
  tab(z, ~ gender + type)
  z <- fillin( z, ~ education + type + gender, xpd = 1.1)
  fit <- lm( income ~ (education+I(education^2)+I(education^3) )* type * gender,
             z, na.action = na.exclude)    # overfitting!
  summary(fit)
  z <- cbind( z, predict(fit, newdata = z, se = TRUE))
  head(z)
  (p <-  xyplot( income ~ education| gender, z, groups = type,
                 subscripts = T,
                 fit = z$fit,
                 lower = z$fit - z$se,
                 upper = z$fit + z$se,
                 layout = c(1,3),
                 auto.key = list(space='right', lines= T, cex = 1.5)))

  p + glayer(panel.fit(...))
  trellis.focus()
  panel.identify(labels= z$occupation)
  trellis.unfocus()
  z$type2 <- with( z, reorder(type,education, mean, na.rm=T))
  gd(3)
  (p <-  xyplot( income ~ education| type2, z, groups = gender,
                 subscripts = T,
                 fit = z$fit,
                 lower = z$fit - z$se,
                 upper = z$fit + z$se,
                 layout = c(1,3),
                 par.strip.text = list(cex = 2),
                 auto.key = list(space='right', lines= T, cex = 1.5)))

  p + glayer( panel.fit(...))
  trellis.focus()
  panel.identify(labels= z$occupation)
  trellis.unfocus()

  ###
  ### With panels^2
  ### need to remove 'col = col.line'
  ###

  z <- Prestige
  z$occ <- rownames(Prestige)
  z$gender <- with(z, cut( women, c(-1,15,50,101),labels = c("Male","Mixed","Female")))
  z$type2 <- with( z, reorder(type,education, mean, na.rm=T))
  tab(z, ~ gender + type2)
  z <- fillin( z, education ~ type + gender, xpd = 1.1)
  fit <- lm( income ~ (education+I(education^2)+I(education^3) )* type * gender,
             z, na.action = na.exclude)    # overfitting!
  summary(fit)
  z <- cbind( z, predict(fit, newdata = z, se = TRUE))
  head(z)
  (p <-  xyplot( income ~ education| gender*type, z,
                 subscripts = T,
                 fit = z$fit,
                 labels = z$occ,
                 lower = z$fit - z$se,
                 upper = z$fit + z$se,
                 auto.key = list(space='right', lines= T, cex = 1.5)))

  p + layer( panel.fit(...))
  p + layer( panel.fit(..., col = 'black', alpha = .1)) + layer(panel.labels(...))
  
  
  

## End(Not run)
## Not run: 
trellis.focus()
panel.identify(labels = rownames(data),rot=-15,col = col.symbol, etc.)
trellis.unfocus()

## End(Not run)

gmonette/spida2 documentation built on Aug. 11, 2024, 7:52 p.m.