Adverse Event Dotplot for R

aedotplot.r

## DEMO FOR SAFETY GRAPHICS

require(lattice)
## The following functions are in the HH package
panel.ae.dotplot <-
 function (x, y, groups, ..., col.AB, pch.AB, lower, upper)
{
    panel.num <- panel.number()
    if (panel.num == 1)
        panel.ae.leftplot(x, y, groups = groups,
            col = col.AB, pch = pch.AB, ...)
    if (panel.num == 2)
        panel.ae.rightplot(x, y, ..., lwd = 6, pch = 16,
                           lower = lower, upper = upper)
}

panel.ae.leftplot <- function(x, y, groups, col.AB, ...)
{
    panel.abline(h = y, lty = 2, lwd = 0, col = gray(0.65))
    panel.superpose(x, y, groups = groups,
                    col = col.AB, ...)
}

panel.ae.rightplot <-
 function(x, y, ..., lwd = 6, lower, upper)
{
    panel.abline(v = 0, lty = 3, lwd = 0, col = 1)
    panel.abline(h = y, lty = 2, lwd = 0, col = gray(0.65))
    panel.segments(lower, y, upper, y, lwd = 2)
    panel.xyplot(x, y, ..., col = 1, cex = 0.7)
    panel.points(lower, y, pch = 3, col = 1, cex = 0.4)
    panel.points(upper, y, pch = 3, col = 1, cex = 0.4)
}


## Version of HH's ae.dotplot that uses risk difference and
## uses proportions instead of percents

riskdiff <- function(data, tx='treat', prefterm='PREF',
                     events='SAE', n='SN',
                     conf.int=.95)
  {
    zcrit <- qnorm((1+conf.int)/2)
    i <- order(data[[prefterm]], data[[tx]])
    data <- data[i,]

    N     <- data[[n]]
    prop  <- data[[events]]/N
    pt    <- data[[prefterm]]
    treat <- data[[tx]]
    if(nrow(data) != 2*length(unique(pt)))
     stop('asymmetric data')
    utrt <- sort(unique(as.character(treat)))

    diff <- prop[treat==utrt[1]] - prop[treat==utrt[2]]
    diff.order <- order(diff)
    pt <- ordered(pt, levels=(pt[treat==utrt[1]])[diff.order])
    diff2 <- as.vector(rbind(diff, diff))
    p1 <- prop[treat==utrt[1]]
    n1 <- N[treat==utrt[1]]
    p2 <- prop[treat==utrt[2]]
    n2 <- N[treat==utrt[2]]

    se <- sqrt(p1*(1-p1)/n1 + p2*(1-p2)/n2)
    se2 <- as.vector(rbind(se, se))
    lower <- diff2 - zcrit*se2
    upper <- diff2 + zcrit*se2
    data.frame(pt, treat, prop, diff=diff2, lower, upper)
}

# Original source of aeanonym: require(HH);
# aeanonym <- read.table(hh("datasets/aedotplot.dat"), header=TRUE, sep=",")

aeanonym <-
structure(list(RAND = structure(c(1L, 2L, 1L, 2L, 1L, 2L, 1L, 
2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 
2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 
2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 
2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L), .Label = c("TREATMENT A", 
"TREATMENT B"), class = "factor"), PREF = structure(c(12L, 12L, 
18L, 18L, 26L, 26L, 33L, 33L, 5L, 5L, 27L, 27L, 6L, 6L, 15L, 
15L, 22L, 22L, 23L, 23L, 31L, 31L, 17L, 17L, 2L, 2L, 3L, 3L, 
13L, 13L, 25L, 25L, 28L, 28L, 14L, 14L, 4L, 4L, 8L, 8L, 19L, 
19L, 21L, 21L, 29L, 29L, 10L, 10L, 20L, 20L, 16L, 16L, 32L, 32L, 
11L, 11L, 1L, 1L, 30L, 30L, 24L, 24L, 9L, 9L, 7L, 7L), .Label = c("ABDOMINAL PAIN", 
"ANOREXIA", "ARTHRALGIA", "BACK PAIN", "BRONCHITIS", "CHEST PAIN", 
"CHRONIC OBSTRUCTIVE AIRWAY", "COUGHING", "DIARRHEA", "DIZZINESS", 
"DYSPEPSIA", "DYSPNEA", "FATIGUE", "FLATULENCE", "GASTROESOPHAGEAL REFLUX", 
"HEADACHE", "HEMATURIA", "HYPERKALEMIA", "INFECTION VIRAL", "INJURY", 
"INSOMNIA", "MELENA", "MYALGIA", "NAUSEA", "PAIN", "RASH", "RESPIRATORY DISORDER", 
"RHINITIS", "SINUSITIS", "UPPER RESP TRACT INFECTION", "URINARY TRACT INFECTION", 
"VOMITING", "WEIGHT DECREASE"), class = "factor"), SAE = c(15L, 
9L, 4L, 9L, 4L, 9L, 2L, 9L, 8L, 11L, 4L, 11L, 9L, 12L, 5L, 12L, 
7L, 12L, 6L, 12L, 6L, 12L, 2L, 14L, 2L, 15L, 1L, 15L, 4L, 16L, 
4L, 17L, 11L, 17L, 6L, 20L, 10L, 23L, 13L, 26L, 12L, 26L, 4L, 
26L, 13L, 28L, 9L, 29L, 12L, 30L, 14L, 36L, 6L, 37L, 8L, 42L, 
20L, 61L, 33L, 68L, 10L, 82L, 23L, 90L, 76L, 95L), SN = c(216L, 
431L, 216L, 431L, 216L, 431L, 216L, 431L, 216L, 431L, 216L, 431L, 
216L, 431L, 216L, 431L, 216L, 431L, 216L, 431L, 216L, 431L, 216L, 
431L, 216L, 431L, 216L, 431L, 216L, 431L, 216L, 431L, 216L, 431L, 
216L, 431L, 216L, 431L, 216L, 431L, 216L, 431L, 216L, 431L, 216L, 
431L, 216L, 431L, 216L, 431L, 216L, 431L, 216L, 431L, 216L, 431L, 
216L, 431L, 216L, 431L, 216L, 431L, 216L, 431L, 216L, 431L)), .Names = c("RAND", 
"PREF", "SAE", "SN"), class = "data.frame", row.names = c(NA, 
-66L))

w <- riskdiff(aeanonym, tx='RAND', prefterm='PREF',
              events='SAE', n='SN')

ae.dotplot <-
  function(data, A.name = 'A', B.name = 'B',
           col.AB = c("red", "blue"), pch.AB = c(16, 17),
           main.title = "Most Frequent On-Therapy Adverse Events Sorted by Risk Difference", 
           main.cex = 1, cex.AB.points = NULL,
           cex.AB.y.scale = 0.6, 
           position.left = c(0, 0, 0.7, 1),
           position.right = c(0.61, 0, 0.98, 1),
           key.y = -0.2, conf.int=0.95) 
{
  r <-
    dotplot(pt ~ prop + diff,
            groups = data$treat,
            data = data, outer = TRUE,
            lower = data$lower, 
            upper = data$upper,
            panel = panel.ae.dotplot, 
            scales = list(x = list(relation = "free",
                          limits = list(range(data$prop), 
                             range(data$lower, data$upper))), 
              y = list(cex = cex.AB.y.scale)),
            A.name = A.name,  B.name = B.name,
            col.AB = col.AB, pch.AB = pch.AB,
            cex.AB.points = cex.AB.points, 
            cex.AB.y.scale = cex.AB.y.scale,
            main = list(main.title, cex = main.cex),
            xlab = list(c(" ", "Less Risk             More Risk"), cex = 0.9), 
            between = list(x = 1), 
            key = list(y = key.y, x = 0.15,
              points = list(col = col.AB, pch = pch.AB),
              text = list(c(A.name, B.name), col = col.AB,
                          cex = 0.9),
              columns = 2, between = 0.5, space = "bottom"))
  r$condlevels[[1]] <- c("Proportion",
                         paste("Risk Difference with",
                          conf.int, "CI"))
  r
}

ae.dotplot(w, A.name='Treatment (N=216)', B.name='Control (N=431)')