Running 'demos.R' [7s/8s]
Running the tests in 'tests/demos.R' failed.
Complete output:
  > library(vcd)
  Loading required package: grid
  > demo(discrete)
  
  
  	demo(discrete)
  	---- ~~~~~~~~
  
  >   #################################################
  >   ## Fitting and Graphing Discrete Distributions ##
  >   #################################################
  > 
  >   data(HorseKicks)
  
  >   barplot(HorseKicks, col = 2,
  +           xlab = "Number of Deaths", ylab = "Number of Corps-Years",
  +           main = "Deaths by Horse Kicks")
  
  >   data(Federalist)
  
  >   barplot(Federalist, col = 2,
  +           xlab = "Occurrences of 'may'", ylab = "Number of Blocks of Text",
  +           main = "'may' in Federalist papers")
  
  >   data(WomenQueue)
  
  >   barplot(WomenQueue, col = 2,
  +           xlab = "Number of women", ylab = "Number of queues",
  +           main = "Women in queues of length 10")
  
  >   data(WeldonDice)
  
  >   barplot(WeldonDice, names = c(names(WeldonDice)[-11], "10+"), col = 2,
  +           xlab = "Number of 5s and 6s", ylab = "Frequency",
  + 	  main = "Weldon's dice data")
  
  >   data(Butterfly)
  
  >   barplot(Butterfly, col = 2,
  +           xlab = "Number of individuals", ylab = "Number of Species",
  +   	main = "Butterfly species im Malaya")
  
  >   ############################
  >   ## Binomial distributions ##
  >   ############################
  > 
  >   par(mfrow = c(1,2))
  
  >   barplot(dbinom(0:10, p = 0.15, size = 10), names = 0:10, col = grey(0.7),
  +           main = "p = 0.15", ylim = c(0,0.35))
  
  >   barplot(dbinom(0:10, p = 0.35, size = 10), names = 0:10, col = grey(0.7),
  +           main = "p = 0.35", ylim = c(0,0.35))
  
  >   par(mfrow = c(1,1))
  
  >   mtext("Binomial distributions", line = 2, cex = 1.5)
  
  >   plot(0:10, dbinom(0:10, p = 0.15, size = 10), type = "b", ylab = "Density",
  +        ylim = c(0, 0.4), main = "Binomial distributions, N = 10", pch = 19)
  
  >   lines(0:10, dbinom(0:10, p = 0.35, size = 10), type = "b", col = 2, pch = 19)
  
  >   lines(0:10, dbinom(0:10, p = 0.55, size = 10), type = "b", col = 4, pch = 19)
  
  >   lines(0:10, dbinom(0:10, p = 0.75, size = 10), type = "b", col = 3, pch = 19)
  
  >   legend(3, 0.4, c("p", "0.15", "0.35", "0.55", "0.75"), lty = rep(1,5), col =
  +          c(0,1,2,4,3), bty = "n")
  
  >   ###########################
  >   ## Poisson distributions ##
  >   ###########################
  > 
  >   par(mfrow = c(1,2))
  
  >   dummy <- barplot(dpois(0:12, 2), names = 0:12, col = grey(0.7), ylim = c(0,0.3),
  +         main = expression(lambda == 2))
  
  >   abline(v = dummy[3], col = 2)
  
  >   diff <- (dummy[3] - dummy[2]) * sqrt(2)/2
  
  >   lines(c(dummy[3] - diff, dummy[3] + diff), c(0.3, 0.3), col = 2)
  
  >   dummy <- barplot(dpois(0:12, 5), names = 0:12, col = grey(0.7), ylim = c(0,0.3),
  +           main = expression(lambda == 5))
  
  >   abline(v = dummy[6], col = 2)
  
  >   diff <- (dummy[6] - dummy[5]) * sqrt(5)/2
  
  >   lines(c(dummy[6] - diff, dummy[6] + diff), c(0.3, 0.3), col = 2)
  
  >   par(mfrow = c(1,1))
  
  >   mtext("Poisson distributions", line = 2, cex = 1.5)
  
  >   #####################################
  >   ## Negative binomial distributions ##
  >   #####################################
  > 
  >   nbplot <- function(p = 0.2, size = 2, ylim = c(0, 0.2))
  +   {
  +     plot(0:20, dnbinom(0:20, p = p, size = size), type = "h", col = grey(0.7),
  +          xlab = "Number of failures (k)", ylab = "Density", ylim = ylim,
  +          yaxs = "i", bty = "L")
  +     nb.mean <- size * (1-p)/p
  +     nb.sd <- sqrt(nb.mean/p)
  +     abline(v = nb.mean, col = 2)
  +     lines(nb.mean + c(-nb.sd, nb.sd), c(0.01, 0.01), col = 2)
  +     legend(14, 0.2, c(paste("p = ", p), paste("n = ", size)), bty = "n")
  +   }
  
  >   par(mfrow = c(3,2))
  
  >   nbplot()
  
  >   nbplot(size = 4)
  
  >   nbplot(p = 0.3)
  
  >   nbplot(p = 0.3, size = 4)
  
  >   nbplot(p = 0.4, size = 2)
  
  >   nbplot(p = 0.4, size = 4)
  
  >   par(mfrow = c(1,1))
  
  >   mtext("Negative binomial distributions for the number of trials to observe n = 2 or n = 4 successes", line = 3)
  
  >   #####################
  >   ## Goodness of fit ##
  >   #####################
  > 
  >   p <- weighted.mean(as.numeric(names(HorseKicks)), HorseKicks)
  
  >   p.hat <- dpois(0:4, p)
  
  >   expected <- sum(HorseKicks) * p.hat
  
  >   chi2 <- sum((HorseKicks - expected)^2/expected)
  
  >   pchisq(chi2, df = 3, lower = FALSE)
  [1] 0.8719411
  
  >   ## or:
  >   HK.fit <- goodfit(HorseKicks)
  
  >   summary(HK.fit)
  
  	 Goodness-of-fit test for poisson distribution
  
                         X^2 df  P(> X^2)
  Likelihood Ratio 0.8682214  3 0.8330891
  
  >   ## Are the dice fair?
  >   p.hyp <- 1/3
  
  >   p.hat <- dbinom(0:12, prob = p.hyp, size = 12)
  
  >   expected <- sum(WeldonDice) * p.hat
  
  >   expected <- c(expected[1:10], sum(expected[11:13]))
  
  >   chi2 <- sum((WeldonDice - expected)^2/expected)
  
  >   G2 <- 2*sum(WeldonDice*log(WeldonDice/expected))
  
  >   pchisq(chi2, df = 10, lower = FALSE)
  [1] 0.0001027878
  
  >   ## Are the data from a binomial distribution?
  >   p <- weighted.mean(as.numeric(names(WeldonDice))/12, WeldonDice)
  
  >   p.hat <- dbinom(0:12, prob = p, size = 12)
  
  >   expected <- sum(WeldonDice) * p.hat
  
  >   expected <- c(expected[1:10], sum(expected[11:13]))
  
  >   chi2 <- sum((WeldonDice - expected)^2/expected)
  
  >   G2 <- 2*sum(WeldonDice*log(WeldonDice/expected))
  
  >   pchisq(chi2, df = 9, lower = FALSE)
  [1] 0.5160827
  
  >   ## or:
  >   WD.fit1 <- goodfit(WeldonDice, type = "binomial", par = list(prob = 1/3, size = 12))
  
  >   WD.fit1$fitted[11] <- sum(predict(WD.fit1, newcount = 10:12))
  
  >   WD.fit2 <- goodfit(WeldonDice, type = "binomial", par = list(size = 12), method = "MinChisq")
  
  >   summary(WD.fit1)
  
  	 Goodness-of-fit test for binomial distribution
  
                        X^2 df     P(> X^2)
  Pearson          37.63908 12 0.0001758723
  Likelihood Ratio 35.10328 10 0.0001198974
  
  >   summary(WD.fit2)
  
  	 Goodness-of-fit test for binomial distribution
  
               X^2 df  P(> X^2)
  Pearson 10.12296 11 0.5193618
  
  >   F.fit1 <- goodfit(Federalist)
  
  >   F.fit2 <- goodfit(Federalist, type = "nbinomial")
  
  >   summary(F.fit1)
  
  	 Goodness-of-fit test for poisson distribution
  
                        X^2 df     P(> X^2)
  Likelihood Ratio 25.24312  5 0.0001250511
  
  >   par(mfrow = c(2,2))
  
  >   plot(F.fit1, scale = "raw", type = "standing")
  
  >   plot(F.fit1, type = "standing")
  
  >   plot(F.fit1)
  
  >   plot(F.fit1, type = "deviation")
  
  >   par(mfrow = c(1,1))
  
  >   plot(F.fit2, type = "deviation")
  
  >   summary(F.fit2)
  
  	 Goodness-of-fit test for nbinomial distribution
  
                        X^2 df  P(> X^2)
  Likelihood Ratio 1.964028  4 0.7423751
  
  >   data(Saxony)
  
  >   S.fit <- goodfit(Saxony, type = "binomial", par = list(size = 12))
  
  >   summary(S.fit)
  
  	 Goodness-of-fit test for binomial distribution
  
                       X^2 df     P(> X^2)
  Likelihood Ratio 97.0065 11 6.978187e-16
  
  >   plot(S.fit)
  
  >   ###############
  >   ## Ord plots ##
  >   ###############
  > 
  >   par(mfrow = c(2,2))
  
  >   Ord_plot(HorseKicks, main = "Death by horse kicks")
  
  >   Ord_plot(Federalist, main = "Instances of 'may' in Federalist papers")
  
  >   Ord_plot(Butterfly, main = "Butterfly species collected in Malaya")
  
  >   Ord_plot(WomenQueue, main = "Women in queues of length 10")
  
  >   par(mfrow = c(1,1))
  
  >   ###############
  >   ## Distplots ##
  >   ###############
  > 
  >   distplot(HorseKicks, type = "poisson")
  
  >   distplot(HorseKicks, type = "poisson", lambda = 0.61)
  
  >   distplot(Federalist, type = "poisson")
  
  >   distplot(Federalist, type = "nbinomial")
  
  >   distplot(Saxony, type = "binomial", size = 12)
  Warning messages:
  1: In summary.goodfit(WD.fit1) :
    Chi-squared approximation may be incorrect
  2: In summary.goodfit(WD.fit2) :
    Chi-squared approximation may be incorrect
  > demo(twoway)
  
  
  	demo(twoway)
  	---- ~~~~~~
  
  >   #####################
  >   ## Fourfold tables ##
  >   #####################
  > 
  >   ### Berkeley Admission Data ###
  >   ###############################
  >   data(UCBAdmissions)
  
  >   ## unstratified
  >   ### no margin is standardized
  >   x <- margin.table(UCBAdmissions, 2:1)
  
  >   fourfold(x, std = "i", extended = FALSE)
  
  >   ### std. for gender
  >   fourfold(x, margin = 1, extended = FALSE)
  
  >   ### std. for both
  >   fourfold(x, extended = FALSE)
  
  >   ## stratified
  >   fourfold(UCBAdmissions, extended = FALSE)
  
  >   fourfold(UCBAdmissions) ## extended plots
  
  >   ## using cotabplot
  >   cotabplot(UCBAdmissions, panel = function(x, condlevels, ...)
  +               fourfold(co_table(x, names(condlevels))[[paste(condlevels, collapse = ".")]],
  +                        newpage = F, return_grob = FALSE, ...)
  +             )
  
  >   ### Coal Miners Lung Data ###
  >   #############################
  >   data(CoalMiners)
  
  >   ## Fourfold display, both margins equated
  >   fourfold(CoalMiners, mfcol = c(3,3))
  
  >   ## Log Odds Ratio Plot
  >   data(CoalMiners, package = "vcd")
  
  >   lor_CM <- loddsratio(CoalMiners)
  
  >   plot(lor_CM)
  
  >   lor_CM_df <- as.data.frame(lor_CM)
  
  >   # fit linear models using WLS
  >   age <- seq(20, 60, by = 5)
  
  >   lmod <- lm(LOR ~ age, weights = 1 / ASE^2, data = lor_CM_df)
  
  >   grid.lines(age, fitted(lmod), gp = gpar(col = "blue"))
  
  >   qmod <- lm(LOR ~ poly(age, 2), weights = 1 / ASE^2, data = lor_CM_df)
  
  >   grid.lines(age, fitted(qmod), gp = gpar(col = "red"))
  
  >   ## Fourfold display, strata equated
  >   fourfold(CoalMiners, std = "ind.max", mfcol = c(3,3))
  
  >   ####################
  >   ## Sieve Diagrams ##
  >   ####################
  > 
  >   ### Hair Eye Color ###
  >   ######################
  >   data(HairEyeColor)
  
  >   ## aggregate over `sex':
  >   (tab <- margin.table(HairEyeColor, 1:2))
         Eye
  Hair    Brown Blue Hazel Green
    Black    68   20    15     5
    Brown   119   84    54    29
    Red      26   17    14    14
    Blond     7   94    10    16
  
  >   ## plot expected values:
  >   sieve(t(tab), sievetype = "expected", shade = TRUE)
  
  >   ## plot sieve diagram:
  >   sieve(t(tab), shade = TRUE)
  
  >   ### Visual Acuity ###
  >   #####################
  >   data(VisualAcuity)
  
  >   attach(VisualAcuity)
  
  >   sieve(Freq ~ right + left,
  +         data = VisualAcuity,
  +         subset = gender == "female",
  +         main = "Unaided distant vision data",
  +         labeling_args = list(set_varnames = c(left = "Left Eye Grade",
  +                                right = "Right Eye Grade")),
  +         shade = TRUE
  +         )
  
  >   detach(VisualAcuity)
  
  >   ### Berkeley Admission ###
  >   ##########################
  > 
  >   ## -> Larger tables: e.g., Cross factors
  >   ### Cross Gender and Admission
  >   data(UCBAdmissions)
  
  >   (tab <- xtabs(Freq ~ Dept + I(Gender : Admit), data = UCBAdmissions))
      I(Gender:Admit)
  Dept Male:Admitted Male:Rejected Female:Admitted Female:Rejected
     A           512           313              89              19
     B           353           207              17               8
     C           120           205             202             391
     D           138           279             131             244
     E            53           138              94             299
     F            22           351              24             317
  
  >   sieve(tab,
  +         labeling_args = list(set_varnames = c("I(Gender:Admit)" = "Gender:Admission",
  +                                Dept = "Department")),
  +         main = "Berkeley Admissions Data",
  +         shade = TRUE
  +         )
  
  >   ## or use extended sieve plots:
  >   sieve(UCBAdmissions, shade = TRUE)
  
  >   ######################
  >   ## Association Plot ##
  >   ######################
  > 
  >   ### Hair Eye Color ###
  >   ######################
  >   data(HairEyeColor)
  
  >   assoc(margin.table(HairEyeColor, 1:2),
  +         labeling_args = list(set_varnames = c(Hair = "Hair Color", Eye = "Eye Color")),
  +         main = "Association Plot")
  
  >   ####################
  >   ## Agreement Plot ##
  >   ####################
  > 
  >   ### Sexual Fun ###
  >   ##################
  >   data(SexualFun)
  
  >   ## Kappa statistics
  >   Kappa(SexualFun)
              value     ASE     z Pr(>|z|)
  Unweighted 0.1293 0.06860 1.885 0.059387
  Weighted   0.2374 0.07832 3.031 0.002437
  
  >   ## Agreement Chart
  >   agreementplot(t(SexualFun), weights = 1)
  
  >   ## Partial Agreement Chart and B-Statistics
  >   (agreementplot(t(SexualFun),
  +                       xlab = "Husband's Rating",
  +                       ylab = "Wife's Rating",
  +                       main = "Husband's and Wife's Sexual Fun")
  +    )
  $Bangdiwala
            [,1]
  [1,] 0.1464624
  
  $Bangdiwala_Weighted
            [,1]
  [1,] 0.4981723
  
  $weights
  [1] 1.0000000 0.8888889
  
  
  >   ### MS Diagnosis data ###
  >   #########################
  >   data(MSPatients)
  
  >   ## use e.g., X11(width = 12), or expand graphics device
  >   agreementplot(t(MSPatients[,,1]), main = "Winnipeg Patients")
  
  >   agreementplot(t(MSPatients[,,2]), main = "New Orleans Patients")
  
  >   ##################
  >   ## Ternary Plot ##
  >   ##################
  > 
  >   ### sample data ###
  >   ###################
  >   (x <- rbind(c(A=10,B=10,C=80),
  +               c(40,30,30),
  +               c(20,60,20)
  +               )
  +    )
        A  B  C
  [1,] 10 10 80
  [2,] 40 30 30
  [3,] 20 60 20
  
  >   ternaryplot(x,
  +               cex = 2,
  +               col = c("black", "blue", "red"),
  +               coordinates = TRUE
  +               )
  
  >   ### Arthritis Treatment Data ###
  >   ################################
  >   data(Arthritis)
  
  >   ## Build table by crossing Treatment and Sex
  >   (tab <- as.table(xtabs(~ I(Sex:Treatment) + Improved, data = Arthritis)))
                  Improved
  I(Sex:Treatment) None Some Marked
    Female:Placebo   19    7      6
    Female:Treated    6    5     16
    Male:Placebo     10    0      1
    Male:Treated      7    2      5
  
  >   ## Mark groups
  >   col <- c("red", "red", "blue", "blue")
  
  >   pch <- c(1, 19, 1, 19)
  
  >   ## plot
  >   ternaryplot(
  +               tab,
  +               col = col,
  +               pch = pch,
  +               cex = 2,
  +               bg = "lightgray",
  +               grid_color = "white",
  +               labels_color = "white",
  +               main = "Arthritits Treatment Data"
  +               )
  
  >   ## legend
  >   grid_legend(0.8, 0.7, pch, col, rownames(tab), title = "GROUP")
  
  >   ### Baseball Hitters Data ###
  >   #############################
  >   data(Hitters)
  
  >   attach(Hitters)
  
  >   colors <- c("black","red","green","blue","red","black","blue")
  
  >   pch <- substr(levels(Positions), 1, 1)
  
  >   ternaryplot(
  +               Hitters[,2:4],
  +               pch = as.character(Positions),
  +               col = colors[as.numeric(Positions)],
  +               main = "Baseball Hitters Data"
  +               )
  
  >   grid_legend(0.8, 0.9, pch, colors, levels(Positions), title = "POSITION(S)")
  
  >   detach(Hitters)
  
  >   ### Lifeboats on the Titanic ###
  >   ################################
  >   data(Lifeboats)
  
  >   attach(Lifeboats)
  The following object is masked from package:datasets:
  
      women
  
  
  >   ternaryplot(
  +               Lifeboats[,4:6],
  +               pch = ifelse(side=="Port", 1, 19),
  +               col = ifelse(side=="Port", "red", "blue"),
  +               id  = ifelse(men/total > 0.1, as.character(boat), NA),
  +               dimnames_position = "edge",
  +               dimnames = c("Men of Crew", "Men passengers", "Women and Children"),
  +               main = "Lifeboats on the Titanic"
  +               )
  
  >   grid_legend(0.8, 0.9, c(1, 19), c("red", "blue"), c("Port", "Starboard"), title = "SIDE")
  
  >   ## Load against time for Port/Starboard boats
  >   plot(launch, total,
  +        pch = ifelse(side == "Port", 1, 19),
  +        col = ifelse(side == "Port", "red", "darkblue"),
  +        xlab = "Launch Time",
  +        ylab = "Total loaded",
  +        main = "Lifeboats on the Titanic"
  +        )
  
  >   legend(as.POSIXct("1912-04-15 01:48:00"), 70,
  +          legend = c("SIDE","Port","Starboard"),
  +          pch = c(NA, 1, 19),
  +          col = c(NA, "red", "darkblue")
  +          )
  
  >   text(as.POSIXct(launch),
  +        total,
  +        labels = as.character(boat),
  +        pos = 3,
  +        offset = 0.3
  +        )
  
  >   abline(lm(total ~ as.POSIXct(launch),
  +             subset = side == "Port"),
  +          col = "red")
  
  >   abline(lm(total ~ as.POSIXct(launch),
  +             subset = side == "Starboard"),
  +          col = "darkblue")
  
  >   detach(Lifeboats)
  > demo(mosaic)
  
  
  	demo(mosaic)
  	---- ~~~~~~
  
  > #####################
  > ## Mosaic Displays ##
  > #####################
  > 
  > #########################
  > ## Hair Eye Color Data ##
  > #########################
  > 
  > data(HairEyeColor)
  
  > ## Basic Mosaic Display ##
  > 
  > HairEye <- margin.table(HairEyeColor, c(1,2))
  
  > mosaic(HairEye, main = "Basic Mosaic Display of Hair Eye Color data")
  
  > ## Hair Mosaic Display with Pearson residuals ##
  > Hair <- margin.table(HairEyeColor,1)
  
  > Hair
  Black Brown   Red Blond 
    108   286    71   127 
  
  > mHair <- as.table(rep(mean(margin.table(HairEyeColor, 1)), 4))
  
  > names(mHair) <- names(Hair)
  
  > mHair
  Black Brown   Red Blond 
    148   148   148   148 
  
  > ## Pearson residuals from Equiprobability model ##
  > 
  > resid <- (Hair - mHair) / sqrt(mHair)
  
  > resid
      Black     Brown       Red     Blond 
  -3.287980 11.343530 -6.329361 -1.726189 
  
  > ## First Step in a Mosaic Display ##
  > 
  > mosaic(Hair, residuals = resid, main = "Hair Color Proportions")
  Error in dimnames(x) <- lapply(d, seq) : 'dimnames' applied to non-array
  Calls: demo ... eval -> eval -> mosaic -> mosaic.default -> strucplot
  Execution halted