#------------------------------------------------------------------------------------- # File name : tas_diamondgraph_figure1.ssc # Date : 12/15/2005 (Version 2) # Authors : Xiuhong Li and Alvaro Muņoz # # Figure 1 in Li X, Buechner J, Tarwater P, Muņoz A. A Diamond-Shaped # Equiponderant Graphical Display of the Effects of Two Categorical Predictors on # Continuous Outcomes. The American Statistician, August 2003, Vol 57, No. 3, 193-199. # # Figure legend:Representations of proportional data (p=0.10,0.25,0.50,and 0.75) as # shaded areas that represent the percentages (10,25,50, and 75) of the # area of each square cell using five alternative methods(A,B,C,D,and E). # # To run code provided, you must have a license to run S-plus. Authors used code to # produce figures in the paper in The American Statistician using S-plus 6.0. # Recipients of code may need to do changes to use in other versions and for # different applications. #------------------------------------------------------------------------------------- #------------------------------------------------------------------------------------- # 1. Obtain the number of categories of the two predictors; # row = number of categories of the predictor on the southwest side of Diamond Graph; # col = number of categories of the predictor on the southeast side of Diamond Graph; #------------------------------------------------------------------------------------- row <- 2; col <- 2; #------------------------------------------------------------------------------------- # 2. Select a square region for plot, and put 6 plots on one page; #------------------------------------------------------------------------------------- par(pty="s",mfrow=c(2,3)); plotfigure1 <- function(panel) # Create a function to plot a panel of figure 1; { #------------------------------------------------------------------------------------- # 3. Assign the coordinate of the lowest point of Diamond Graph to (0,0), and define # the ranges for x and y axes based on the number of categories of two predictors. #------------------------------------------------------------------------------------- plot(0,0, xlim=c(0-row, col), ylim=c(0,row+col),axes=F, type="n", xlab="", ylab=""); #------------------------------------------------------------------------------------- # 4. Obtain the coordinate of the center of each square cell (cx,cy) with diagonal= 2; # i.e., the distance from the center to a vertex of a square cell equals to 1. #------------------------------------------------------------------------------------- for (i in 1:row) { if (i == 1){ cx <- seq(0,(col-1) ) cy <- seq(1,col) } else{ cx <- c(cx, seq( (1-i),(col-i) ) ) cy <- c(cy, seq(i,(col+i-1)) ) } } #------------------------------------------------------------------------------------- # 5. Read in the outcomes, which start from the value in the cell 1 (at bottom) # to cell 2, to cell 3, and finally to cell 4 (at top), as depicted below in # the 2x2 graph. # /4\ <- top # / \/\ # \3/\2/ # \1/ <- bottom # # Empty cells (with no subjects) are assigned -1 and will not be displayed. # Only positive outcomes will be plotted as polygons in square cells. #------------------------------------------------------------------------------------- outcome <- c(0.1,0.25,0.5,0.75); #------------------------------------------------------------------------------------- # 6.If data are not proportional, convert the outcomes to proportions; # to get relative risks,each value is divided by the highest one within that outcome; # to get excess risks, the difference between each value and the minimum outcome is # divided by the difference between the maximum and the minimum outcome. # If data are proportional, no need to convert, which is the case in figure 1. #------------------------------------------------------------------------------------- p <- outcome #------------------------------------------------------------------------------------- # 7. Calculate the 'scale' based on the number of categories of two predictors; # We have found scaling all proportions (p) by 'scale' useful in the # graphical illustration when p is high. # If a page has only one graph, we use scale = 1.005 - 0.005(row+col). # But figure 1 has 6 panels; the size of each square in each panel is about # the size of the square of a 3x6 plot, so we used # scale =1.005 - 0.005*(3+6)=0.96 for A-D, and scale=0.99 for E. #------------------------------------------------------------------------------------- scale <- 1.005 - 0.005*(3+6); if (panel=="E") scale <- 0.99; p[p<1] <- p[p<1] *scale p[p==1] <- p[p==1]*scale+0.005 #------------------------------------------------------------------------------------- # 8. Draw the outline of the square cells. # 1) draw 45 degree lines from southwest to northeast # 2) draw -45 degree lines from southeast to northwest #------------------------------------------------------------------------------------- for (i in 0:row) lines(c(-i,col-i),c(i,col+i),lwd=0.01) for (i in 0:col) lines(c(i, -row+i), c(i, row+i),lwd=0.01) #------------------------------------------------------------------------------------- # 9. Calculate the coordinates of the vertices of inside polygons. The polygon below # was an example of using method C; start clockwise from northwest vertex as (x1, y1). # x1,y1 _____ x2,y2 # / \ # x6,y6 / \ x3,y3 # \ / # \_____/ # x5,y5 x4,y4 #------------------------------------------------------------------------------------- if (panel=="A") { x1<-cx x2<-cx + sqrt(p) x3<-cx x4<-cx - sqrt(p) y1<-cy + sqrt(p) y2<-cy y3<-cy - sqrt(p) y4<-cy } else if (panel=="B"){ x1<-cx - sqrt(1-p) x2<-cx + sqrt(1-p) x3<-cx + 1 x4<-cx + sqrt(1-p) x5<-cx - sqrt(1-p) x6<-cx - 1 y1<-cy + (1-sqrt(1-p)) y2<-cy + (1-sqrt(1-p)) y3<-cy y4<-cy - (1-sqrt(1-p)) y5<-cy - (1-sqrt(1-p)) y6<-cy } else if (panel=="C") { x1<- cx - (1-p)/2 x2<- cx + (1-p)/2 x3<- cx + (1+p)/2 x4<- cx + (1-p)/2 x5<- cx - (1-p)/2 x6<- cx - (1+p)/2 y1<- cy + p y2<- cy + p y3<- cy y4<- cy - p y5<- cy - p y6<- cy } else if (panel=="D"){ x1<- cx x2<- cx + 1 x3<- cx x4<- cx - 1 y1<- cy + p y2<- cy y3<- cy - p y4<- cy } else if (panel=="E"){ x1<- cx - (1-p) x2<- cx + (1-p) x3<- cx + p x4<- cx + (1-p) x5<- cx - (1-p) x6<- cx - p y1<- cy + p y2<- cy + p y3<- cy y4<- cy - p y5<- cy - p y6<- cy } #------------------------------------------------------------------------------------- # 10. Plot a polygon inside each square cell, and then print the value of the outcome # at the center in each square cell. # # Here we use cex= 0.7*scale for the labels. Other applications may # need to use a different factor so that labels can be read easily. #------------------------------------------------------------------------------------- n <- row*col for (i in 1:n) { if (panel=="A" | panel=="D") polygon(c(x1[i],x2[i],x3[i],x4[i]), c(y1[i],y2[i],y3[i],y4[i]), density=-.1, border=F, col=5) else if (panel=="B" | panel=="C" | panel=="E") polygon(c(x1[i],x2[i],x3[i],x4[i],x5[i],x6[i]), c(y1[i],y2[i],y3[i],y4[i],y5[i],y6[i]), density=-.1, border=F, col=5) if (panel==""){ # In Window system, S-plus needs 0.03 to be added to cy so that lable is properly centered text(cx[i], cy[i]+0.03, outcome[i]*100, cex=0.7*scale, lwd=2); } } title(main=panel) }# end of function plotfigure1; # Call the function to plot 6 panels in figure 1 plotfigure1("") plotfigure1("A") plotfigure1("B") plotfigure1("C") plotfigure1("D") plotfigure1("E")