rm(list=ls(all=TRUE))
#you will need the following libraries:
library(arm)
library(boot)
library(plotrix) #for calculating standard error of mean

#1) Case-level analysis, dropping cases from before 1961

#case-level analysis
#keep only cases from 1961
case.data <- read.dta("songer_caselevel_1925_2002_vote_data.dta", convert.underscore=T)
case.data <- case.data[case.data$year >=1961,]
attach.all(case.data, overwrite=T)
case <- casenum
year <- year
circuit <- ifelse(circuit == 0, 12, circuit)
issue <- geniss
pvote <- panel.vote
lower.vote <- lower.vote
j1.dem <- ifelse(j1presparty == 100,1,0)
j2.dem <- ifelse(j2presparty == 100,1,0)
j3.dem <- ifelse(j3presparty == 100,1,0)
j1.ideal <-j1jcs
j2.ideal <-j2jcs
j3.ideal <-j3jcs
j1.code <- j1code
j2.code <- j2code
j3.code <- j3code
mixed <- mixed

keep.data <- data.frame(case, year, circuit, issue, pvote, j1.dem, j2.dem, j3.dem,
    j1.ideal, j2.ideal, j3.ideal, lower.vote, j1.code, j2.code, j3.code)
detach()
attach.all(keep.data, overwrite = T)
keep.data <- keep.data[order(year,circuit, case),]
detach()
attach.all(keep.data)
keep.data$circuit.year <- paste(circuit,year)
detach()
attach.all(keep.data)
unique.circuit.year <- unique(circuit.year)
unique.year <- unique(year)
loop <- length(unique.circuit.year)

#attach Songer weights to data

#get weights to estimate overall conservative rate; create table like this:

#                sample of circuit          universe of circuit      
#circuit     # cases     proportion   # cases   proportion   weight   
#_________________________________________________________________    
#01              15      .1              095     .049       0.49     
#02              15      .1              329     .170       1.70     
#03              15      .1              116     .060       0.60     
#04              15      .1              099     .051       0.51     
#05              15      .1              175     .091       0.91     
#06              15      .1              222     .115       1.15     
#07              15      .1              081     .042       0.42     
#08              15      .1              330     .171       1.71     
#09              15      .1              289     .150       1.50     
#DC              15      .1              196     .101       1.01     
#_______________________________________________________________      
#total          150     1.0              1932    1.0      

#use circuit-year as unit of analysis
songer.weights <- read.dta("songer_weight_data.dta", convert.underscore=T)
songer.weights <- songer.weights[songer.weights$year >= 1961,]
songer.weights$circuit <- ifelse(songer.weights$circuit == 0, 12, songer.weights$circuit)
detach()
attach.all(songer.weights, overwrite=T)

songer.weights <- songer.weights[order(year,circuit),]
songer.weights$circuit.year <- paste(songer.weights$circuit,songer.weights$year)#for removing cases below
detach()
attach.all(songer.weights, overwrite=T)
n.cases.circuit.universe <- cases #number of cases in universe, by circuit-year

#nr <- nrow(songer.weights)#shortcut for loops below
#since keep.data has fewer circuit-years than songer, use length of keep.data for looping
n.cases.circuit.sample <- rep(NA, loop)#number of cases in a circuit, in a given year
n.cases.year.sample <- rep(NA, loop)#number of cases total, in a given year
sample.proportion <- rep(NA, loop)#propotion of all in-sample cases, for a given year
universe.proportion <- rep(NA, loop)#proportion of all cases in universe, for a given year
#create condensed dataset to merge workload onto circuit-year level data below.
caseload.circuit.year <- data.frame(cbind(year, circuit, cases))


#SAMPLE MEASURES:  get number of cases and proportions from actual Songer data, by circuit and year
detach()
attach.all(keep.data)

#get number of cases, per circuit-year for sample 
for (i in 1:loop){ #loop over circuit year
     n.cases.circuit.sample[i] <- length(case[circuit.year == unique.circuit.year[i]])
            }
            
#first, create vector based on unique-year, then we will attach to circuit-year data.
n.cases.byyear.sample <- rep(NA, length(unique.year))
for(i in 1:length(unique.year)){
    n.cases.byyear.sample[i] <-  sum(length(case[year == unique.year[i]]))
    }
n.cases.byyear.data <- data.frame(n.cases.byyear.sample, unique.year)#for merging to circuit-year level
n.cases.byyear.data$year <- n.cases.byyear.data$unique.year
n.cases.byyear.data <- n.cases.byyear.data[,c(1,3)]
    
#expand to circuit-year    
detach()
attach.all(songer.weights, overwrite=T)
songer.weights <- data.frame(cbind(songer.weights, n.cases.circuit.sample))
songer.weights <- merge(songer.weights, n.cases.byyear.data)
#get sample proportion
songer.weights$sample.proportion  <- songer.weights$n.cases.circuit.sample/songer.weights$n.cases.byyear.sample
detach()
attach.all(songer.weights, overwrite=T)

###############################
##UNIVERSE MEASURES

#get total number of cases in each year:
#first, create vector based on unique-year, then we will attach to circuit-year data.
n.cases.byyear.universe <- rep(NA, length(unique.year))
for (i in 1:length(n.cases.byyear.universe)){
    n.cases.byyear.universe[i] <- sum(cases[year==unique.year[i]])
    }
#expand to circuit-year    
universe.year.data <- data.frame(n.cases.byyear.universe, unique.year)#for merging
universe.year.data$year <- universe.year.data$unique.year
universe.year.data <- universe.year.data[,c(1,3)]
songer.weights <- merge(songer.weights, universe.year.data) 
 #get universe proportion
songer.weights$universe.proportion  <- songer.weights$cases/songer.weights$n.cases.byyear.universe
#GET WEIGHTS = DIVIDE UNIVERSE PROPORTION BY SAMPLE PROPORTION
songer.weights$weight <- songer.weights$universe.proportion / songer.weights$sample.proportion
detach()
attach.all(songer.weights, overwrite=T)

#merge weights with case-level data

keep.data <- merge(keep.data, songer.weights)
detach()
attach.all(keep.data)
keep.data <- keep.data[order(year,circuit, case),]
detach()
attach.all(keep.data)

#set up dataframe to merge weights with judge-level data below
weights.to.merge <- data.frame(cbind(year, circuit,weight))
weights.to.merge <- unique(weights.to.merge)#every obs within a circuit-year has same weight,so collapse dataset into circuit-year level

################################################
#Descriptive Analysis of Panel Effects:  All Years

#code party configurations of panel
keep.data$ddd <- ifelse(j1.dem == 1 & j2.dem ==1 & j3.dem ==1, 1, 0)
keep.data$rrr   <- ifelse(j1.dem == 0 & j2.dem == 0 & j3.dem == 0, 1, 0)
keep.data$ddr <- ifelse(j1.dem + j2.dem + j3.dem == 2, 1,0)
keep.data$rrd <-  ifelse(j1.dem + j2.dem + j3.dem == 1, 1,0)

#pool within 5-year intervals, except at end

keep.data$interval <-  cut(year, c(1960, 1965,1970,1975,1980,1985,1990, 1995, 1999, 2003),
  labels=c("1961-1965","1966-1970", "1971-1975","1976-1980","1981-1985","1986-1990", "1991-1995", "1996-1999",
    "2000-2002"))
#create codes for different case types
keep.data$criminal <- ifelse(issue == 1,1,0)
keep.data$civil <- ifelse(issue >=2 & issue <=5,1,0)#combine civil rights, First Amend, due process and privacy
keep.data$econ <- ifelse(issue == 6 | issue == 7,1,0)#combine labor, econ
detach()
attach.all(keep.data)    
    
x.label <- c("1961-65","1966-70", "1971-75","1976-80","1981-85","1986-90", "1991-95", "1996-99",
    "2000-02")

unique.interval <- unique(interval)
loop <- length(unique.interval)

#get number of cases in each interval
n.cases <- rep(NA, loop)
for (i in 1:loop){
    keep <- interval == unique.interval[i]
    n.cases[i] <- length(pvote[keep])
    }
    
rrr.all.n <- rep(NA, loop)
rrr.all.rate <- rep(NA, loop)
rrr.all.ci <- matrix(nrow = loop, ncol =2)
ddd.all.n <- rep(NA, loop)
ddd.all.rate <- rep(NA, loop)
ddd.all.ci <- matrix(nrow = loop, ncol =2)
rrd.all.n <- rep(NA, loop)
rrd.all.rate <- rep(NA, loop)
rrd.all.ci <- matrix(nrow = loop, ncol =2)
ddr.all.n <- rep(NA, loop)
ddr.all.rate <- rep(NA, loop)
ddr.all.ci <- matrix(nrow = loop, ncol =2)
rrr.all.test <- rep(NA, loop)


rrr.all.mixed.rate <- rep(NA, loop)
rrr.all.mixed.rate.ci <- matrix(nrow = loop, ncol =2)
rrr.all.mixed.rate <- rep(NA, loop)
rrr.all.mixed.rate.ci <- matrix(nrow = loop, ncol =2)


#for bootstrapping
wmean <- function(x,k) { #x is cbind(pvote,weight)
  weighted.mean(x[k,1],x[k,2])
}

#get results
for (i in 1:loop){
  keep <- interval == unique.interval[i] #use this for intervals
  
#RRR
        keep2 <- rrr==1
    rrr.all.n[i] <- length(pvote[keep & keep2])
    rrr.all.rate[i] <-  weighted.mean(pvote[keep & keep2], weight[keep & keep2])
    rrr.all.boot <-  boot(cbind(pvote[keep & keep2], weight[keep & keep2]), wmean, R = 1000, stype = "i", sim = "ordinary")
    rrr.all.ci[i,] <- quantile(rrr.all.boot$t,probs=c(.025,.975))
    rrr.all.test[i] <- mean(pvote[keep & keep2])
   
###############################
#DDD
          keep2 <- ddd==1
    ddd.all.n[i] <- length(pvote[keep & keep2])
    ddd.all.rate[i] <-  weighted.mean(pvote[keep & keep2], weight[keep & keep2])
    ddd.all.boot <-  boot(cbind(pvote[keep & keep2], weight[keep & keep2]), wmean, R = 1000, stype = "i", sim = "ordinary")
    ddd.all.ci[i,] <- quantile(ddd.all.boot$t,probs=c(.025,.975))

#########################################
#RRD
          keep2 <- rrd==1
    rrd.all.n[i] <- length(pvote[keep & keep2])
    rrd.all.rate[i] <-  weighted.mean(pvote[keep & keep2], weight[keep & keep2])
    rrd.all.boot <-  boot(cbind(pvote[keep & keep2], weight[keep & keep2]), wmean, R = 1000, stype = "i", sim = "ordinary")
    rrd.all.ci[i,] <- quantile(rrd.all.boot$t,probs=c(.025,.975))
 
#########################################
#DDR
        keep2 <- ddr==1
    ddr.all.n[i] <- length(pvote[keep & keep2])
    ddr.all.rate[i] <-  weighted.mean(pvote[keep & keep2], weight[keep & keep2])
    ddr.all.boot <-  boot(cbind(pvote[keep & keep2], weight[keep & keep2]), wmean, R = 1000, stype = "i", sim = "ordinary")
    ddr.all.ci[i,] <- quantile(ddr.all.boot$t,probs=c(.025,.975))
  }

#plotting shortcuts
x.axis <- c(1:length(unique.interval))
axis.text <- 1.6
rrr.color <- "red"
ddd.color <- "blue"
rrd.color <- "purple"
ddr.color <- "darkcyan"

#for each graph: bracket intervals using dotted vertical lines
  #then stagger one point and the other so can present confidence intervals
point.stagger <- .1  
right.offset <- x.axis + .5 + point.stagger #add .5 so points are centered in each interval
left.offset <- x.axis + .5 -  point.stagger
point.size <- .9
axis.text <- 1.2
x.tick.marks <- c(1:10)
x.lim <- c(1,10)
y.lim <- c(.15, .7)



##########################
#CREATE GRAPHS

#1)Figure 4: DO  RRR vs DDD Vote differently over time? is it constant

pdf("figure4.pdf", height = 3.5, width =10)    

par(mar = c(3,3.5,1,.5))
plot(right.offset, rrr.all.rate, type = "n", xlim = x.lim, ylim = y.lim, main ="",
    xlab = "", ylab = "", axes = F, xaxs = "i", yaxs = "i", col = rrr.color)
segments(right.offset, rrr.all.ci[,1], right.offset, rrr.all.ci[,2], col = rrr.color, lty = 2)#add 95% confidence intervals
segments(left.offset, ddd.all.ci[,1], left.offset, ddd.all.ci[,2], col = ddd.color)
points(left.offset, ddd.all.rate, type = "p", col = ddd.color,  pch = 19, cex = point.size)#add estimates
points(right.offset, rrr.all.rate, type = "p", col = rrr.color,  pch = 21, bg = "white", cex = point.size)#bg plots points over lines
abline(v=x.axis, lty = 2, lwd = .5, col = "light gray")
axis(1, at = x.tick.marks, label = NA, tick = T)#add tick marks at interval breaks
axis(1, at = x.axis + .5, label = x.label, tick = F, mgp = c(2,.7,0),cex.axis = axis.text)#add labels in middle of intervals
axis(2, at = seq(0,1, by = .1), label = seq(0,100, by = 10), las =2,
    cex.axis = axis.text, mgp = c(2,.7,0))
mtext("Percent Liberal", 2, line = 2.5, cex = axis.text)
#add legend
top.point.height <- .68
bottom.point.height <- top.point.height - .05
x.location <- 6.1
points(x.location, top.point.height, type = "p", col = ddd.color,  pch = 19, cex = point.size)
text(x.location + .05, top.point.height, "DDD panels",col = ddd.color, adj = 0) 
points(x.location + .1,bottom.point.height, type = "p", col = rrr.color, cex = point.size, pch = 21, bg = "white")
text(x.location + .15, bottom.point.height, "RRR panels",col = rrr.color, adj = 0) 

dev.off()

##############################

#Figure 5: How different are unified panels from mixed panels?
  #plot both DDD vs. DDR and RRR vs. RRD


pdf("figure5.pdf", height = 7, width =10)    

#DDD vs. ddr
par(mfrow = c(2,1), mar = c(3,3.5,2.5,.5))
plot(right.offset, ddr.all.rate, type = "n", xlim = x.lim, ylim = y.lim, main ="",
    xlab = "", ylab = "", axes = F, xaxs = "i", yaxs = "i", col = ddr.color)
segments(right.offset, ddr.all.ci[,1], right.offset, ddr.all.ci[,2], col = ddr.color, lty = 2)#add 95% confidence intervals
segments(left.offset, ddd.all.ci[,1], left.offset, ddd.all.ci[,2], col = ddd.color)
points(left.offset, ddd.all.rate, type = "p", col = ddd.color,  pch = 19, cex = point.size)#add estimates
points(right.offset, ddr.all.rate, type = "p", col = ddr.color,  pch = 21, bg = "white", cex = point.size)#bg plots points over lines
abline(v=x.axis, lty = 2, lwd = .5, col = "light gray")
axis(1, at = x.tick.marks, label = NA, tick = T)#add tick marks at interval breaks
axis(1, at = x.axis + .5, label = x.label, tick = F, mgp = c(2,.7,0),cex.axis = axis.text)#add labels in middle of intervals
axis(2, at = seq(0,1, by = .1), label = seq(0,100, by = 10), las =2,
    cex.axis = axis.text, mgp = c(2,.7,0))
mtext("Percent Liberal", 2, line = 2.5, cex = axis.text)
mtext("Majority Democratic Panels", 3, line = 1, cex = 1.4)
#add legend
top.point.height <- .68
bottom.point.height <- top.point.height - .05
x.location <- 6.1
points(x.location, top.point.height, type = "p", col = ddd.color,  pch = 19, cex = point.size)
text(x.location + .05, top.point.height, "DDD panels",col = ddd.color, adj = 0) 
points(x.location + .1,bottom.point.height, type = "p", col = ddr.color, cex = point.size, pch = 21, bg = "white")
text(x.location + .15, bottom.point.height, "DDR panels",col = ddr.color, adj = 0) 

####################################
#2)#RRR vs. RRD

plot(right.offset, rrr.all.rate, type = "n", xlim = x.lim, ylim = y.lim, main ="",
    xlab = "", ylab = "", axes = F, xaxs = "i", yaxs = "i", col = rrr.color)
segments(right.offset, rrr.all.ci[,1], right.offset, rrr.all.ci[,2], col = rrr.color, lty = 2)#add 95% confidence intervals
segments(left.offset, rrd.all.ci[,1], left.offset, rrd.all.ci[,2], col = rrd.color)
points(left.offset, rrd.all.rate, type = "p", col = rrd.color,  pch = 19, cex = point.size)#add estimates
points(right.offset, rrr.all.rate, type = "p", col = rrr.color,  pch = 21, bg = "white", cex = point.size)#bg plots points over lines
abline(v=x.axis, lty = 2, lwd = .5, col = "light gray")
axis(1, at = x.tick.marks, label = NA, tick = T)#add tick marks at interval breaks
axis(1, at = x.axis + .5, label = x.label, tick = F, mgp = c(2,.7,0),cex.axis = axis.text)#add labels in middle of intervals
axis(2, at = seq(0,1, by = .1), label = seq(0,100, by = 10), las =2,
    cex.axis = axis.text, mgp = c(2,.7,0))

mtext("Majority Republican Panels", 3, line = 1, cex = 1.4)   
mtext("Percent Liberal", 2, line = 2.5, cex = axis.text)
#add legend
top.point.height <- .68
bottom.point.height <- top.point.height - .05
x.location <- 6.1
points(x.location, top.point.height, type = "p", col = rrd.color,  pch = 19, cex = point.size)
text(x.location + .05, top.point.height, "RRD panels",col = rrd.color, adj = 0) 
points(x.location + .1,bottom.point.height, type = "p", col = rrr.color, cex = point.size, pch = 21, bg = "white")
text(x.location + .15, bottom.point.height, "RRR panels",col = rrr.color, adj = 0) 

dev.off()


########################################
####################################
##JUDGE-LEVEL Data

judge.data <- read.dta("songer_judgelevel_1925_2002_vote_data.dta", convert.underscore=T)
detach()
judge.data <-judge.data[judge.data$year>=1961,]#keep only cases decided after 1960
attach.all(judge.data, overwrite=T)
#note: for this data, j1 is the unit of analysis; use j2 and j3 to measure patterns in voting.
issue <- geniss
circuit <- ifelse(circuit ==0, 12, circuit)#recode DC
circuit.year <- paste(circuit, year)
year <- year
j1.dem <- ifelse(j1presparty == 100,1,0)
j1.code <- j1code
j1.ideal <- j1jcs
j1.vote <- j1vote
j2.dem <- ifelse(j2presparty == 100,1,0)
j2.code <- j2code
j2.ideal <- j2jcs
j2.vote <- ifelse(j2vote==1,0,1)
j3.dem <- ifelse(j3presparty == 100,1,0)
j3.code <- j3code
j3.ideal <- j3jcs
j3.vote <- ifelse(j3vote==1,0,1)
lower.vote <- ifelse(lower.vote,0,1)

keep.data.judge <- data.frame(case = judge.data$casenum, year, circuit, circuit.year,
    issue, j1.dem, j2.dem, j3.dem, j1.ideal, j2.ideal, j3.ideal, j1.code, j2.code, 
    j3.code, j1.vote, j2.vote, j3.vote, lower.vote)
detach()    
attach.all(keep.data.judge, overwrite=T)
keep.data.judge <- keep.data.judge[order(year,circuit, case),]    
detach()    
attach.all(keep.data.judge, overwrite=T)


#attach weights to data
keep.data.judge <- merge(keep.data.judge, weights.to.merge)
detach()
attach.all(keep.data.judge)

#code panels based on j2 an j3
keep.data.judge$rr.panel <- ifelse(j2.dem + j3.dem == 0,1,0)
keep.data.judge$dd.panel <- ifelse(j2.dem + j3.dem == 2,1,0)
keep.data.judge$rd.panel <- ifelse(j2.dem + j3.dem == 1,1,0)


#pool within 5-year intervals, except at end

keep.data.judge$interval <-  cut(year, c(1960, 1965,1970,1975,1980,1985,1990, 1995, 1999, 2003),
  labels=c("1961-1965","1966-1970", "1971-1975","1976-1980","1981-1985","1986-1990", "1991-1995", "1996-1999",
    "2000-2002"))
#create codes for different case types
keep.data.judge$criminal <- ifelse(issue == 1,1,0)
keep.data.judge$civil <- ifelse(issue >=2 & issue <=5,1,0)#combine civil rights, First Amend, due process and privacy
keep.data.judge$econ <- ifelse(issue == 6 | issue == 7,1,0)#combine labor, econ
detach()
attach.all(keep.data.judge)    

x.label <- c("1961-65","1966-70", "1971-75","1976-80","1981-85","1986-90", "1991-95", "1996-99",
    "2000-02")

unique.interval <- unique(interval)
loop <- length(unique.interval)


dem.dd.rate <- rep(NA, loop)  #demjudges with dd colleagues
dem.dd.ci <- matrix(nrow = loop, ncol =2)
dem.rd.rate <- rep(NA, loop)  #dem judges with dr colleagues
dem.rd.ci <- matrix(nrow = loop, ncol =2)
dem.rr.rate <- rep(NA, loop) #dem judges with rr colleagues
dem.rr.ci <- matrix(nrow = loop, ncol =2)
gop.rr.rate  <- rep(NA, loop)#gop judges with rr colleagues
gop.rr.ci  <- matrix(nrow = loop, ncol =2)
gop.rd.rate <- rep(NA, loop) #gop judges with rd colleagues
gop.rd.ci <- matrix(nrow = loop, ncol =2)
gop.dd.rate <- rep(NA, loop) #gop judges with dd colleagues
gop.dd.ci <- matrix(nrow = loop, ncol =2)

#get results
for (i in 1:loop){
   keep1 <- interval == unique.interval[i] #use this for intervals
        
        keep2 <- j1.dem == 1 & dd.panel == 1
    dem.dd.rate[i] <-  weighted.mean(j1.vote[keep1 & keep2], weight[keep1 & keep2])
    dem.dd.boot <-  boot(cbind(j1.vote[keep1 & keep2], weight[keep1 & keep2]), wmean, R = 1000, stype = "i", sim = "ordinary")
    dem.dd.ci[i,] <-  quantile(dem.dd.boot$t,probs=c(.025,.975))
       
        keep2 <- j1.dem == 1 & rd.panel == 1
    dem.rd.rate[i] <-  weighted.mean(j1.vote[keep1 & keep2], weight[keep1 & keep2])
    dem.rd.boot <-  boot(cbind(j1.vote[keep1 & keep2], weight[keep1 & keep2]), wmean, R = 1000, stype = "i", sim = "ordinary")
    dem.rd.ci[i,] <-  quantile(dem.rd.boot$t,probs=c(.025,.975))
      
        keep2 <- j1.dem == 1 & rr.panel == 1
    dem.rr.rate[i] <-  weighted.mean(j1.vote[keep1 & keep2], weight[keep1 & keep2])
    dem.rr.boot <-   boot(cbind(j1.vote[keep1 & keep2], weight[keep1 & keep2]), wmean, R = 1000, stype = "i", sim = "ordinary")
    dem.rr.ci[i,] <-  quantile(dem.rr.boot$t,probs=c(.025,.975)) 
               
        keep2 <- j1.dem == 0 & rr.panel == 1
    gop.rr.rate[i] <-  weighted.mean(j1.vote[keep1 & keep2], weight[keep1 & keep2])
    gop.rr.boot <-  boot(cbind(j1.vote[keep1 & keep2], weight[keep1 & keep2]), wmean, R = 1000, stype = "i", sim = "ordinary")
    gop.rr.ci[i,] <-  quantile(gop.rr.boot$t,probs=c(.025,.975))
   
        keep2 <- j1.dem == 0 & rd.panel == 1
    gop.rd.rate[i] <-  weighted.mean(j1.vote[keep1 & keep2], weight[keep1 & keep2])
    gop.rd.boot <-  boot(cbind(j1.vote[keep1 & keep2], weight[keep1 & keep2]), wmean, R = 1000, stype = "i", sim = "ordinary")
    gop.rd.ci[i,] <-  quantile(gop.rd.boot$t,probs=c(.025,.975))
   
        keep2 <- j1.dem == 0 & dd.panel == 1
    gop.dd.rate[i] <-  weighted.mean(j1.vote[keep1 & keep2], weight[keep1 & keep2])
    gop.dd.boot <-  boot(cbind(j1.vote[keep1 & keep2], weight[keep1 & keep2]), wmean, R = 1000, stype = "i", sim = "ordinary")
    gop.dd.ci[i,] <-  quantile(gop.dd.boot$t,probs=c(.025,.975))
}
############
#Figure 6: 

#plot shortcuts
rr.color <- "red"
rd.color <- "dark green"
dd.color <- "blue"


pdf("figure6.pdf", height = 7, width =10)    
#1) DEMOCRATS
par(mfrow = c(2,1), mar = c(4,3.5,2,.5))
plot(right.offset, dem.dd.rate, type = "n", xlim = x.lim, ylim = y.lim, main ="",
    xlab = "", ylab = "", axes = F, xaxs = "i", yaxs = "i", col = ddr.color)
segments(right.offset,  dem.rr.ci[,1], right.offset, dem.rr.ci[,2], col = rr.color)
segments(x.axis +.5, dem.rd.ci[,1], x.axis +.5, dem.rd.ci[,2], col = rd.color)#add 95% confidence intervals
segments(left.offset,  dem.dd.ci[,1], left.offset, dem.dd.ci[,2], col = dd.color)#add 95% confidence intervals

points(right.offset, dem.rr.rate, type = "p", col = rr.color,  pch = 21, bg = "white", cex = point.size)#bg plots points over lines
points(x.axis+.5, dem.rd.rate, type = "p", col = rd.color,  pch = 23, bg = "white", cex = point.size)#add estimates
points(left.offset, dem.dd.rate, type = "p", col = dd.color,  pch = 19, cex = point.size)#add estimates

abline(v=x.axis, lty = 2, lwd = .5, col = "light gray")
axis(1, at = x.tick.marks, label = NA, tick = T)#add tick marks at interval breaks
axis(1, at = x.axis + .5, label = x.label, tick = F, mgp = c(2,.7,0),cex.axis = axis.text)#add labels in middle of intervals
axis(2, at = seq(0,1, by = .1), label = seq(0,100, by = 10), las =2,
    cex.axis = axis.text, mgp = c(2,.7,0))
mtext("Percent Liberal", 2, line = 2.5, cex = axis.text)
mtext("Democratic Judges", 3, line = .2, cex = 1.4)
#add legend
top.point.height <- .68
middle.point.height <- top.point.height - .05
bottom.point.height <- top.point.height - .1
x.location <- 6.5
points(x.location, top.point.height, type = "p", col = dd.color,  pch = 19, cex = point.size)
text(x.location + .05, top.point.height, "DD",col = dd.color, adj = 0) 
points(x.location + .1,middle.point.height, type = "p", col = rd.color, cex = point.size, pch = 23, bg = "white")
text(x.location + .15, middle.point.height, "RD",col = rd.color, adj = 0) 
points(x.location + .2,bottom.point.height, type = "p", col = rr.color, cex = point.size, pch = 21)
text(x.location + .25, bottom.point.height, "RR",col = rr.color, adj = 0) 

#2) REPUBLICAN JUDGES

plot(right.offset, gop.dd.rate, type = "n",  x.lim, ylim = y.lim, main ="",
    xlab = "", ylab = "", axes = F, xaxs = "i", yaxs = "i", col = ddr.color)
segments(right.offset,  gop.rr.ci[,1], right.offset, gop.rr.ci[,2], col = rr.color)
segments(x.axis +.5, gop.rd.ci[,1], x.axis +.5, gop.rd.ci[,2], col = rd.color)#add 95% confidence intervals
segments(left.offset,  gop.dd.ci[,1], left.offset, gop.dd.ci[,2], col = dd.color)#add 95% confidence intervals

points(right.offset, gop.rr.rate, type = "p", col = rr.color,  pch = 21, bg = "white", cex = point.size)#bg plots points over lines
points(x.axis+.5, gop.rd.rate, type = "p", col = rd.color,  pch = 23, bg = "white", cex = point.size)#add estimates
points(left.offset, gop.dd.rate, type = "p", col = dd.color,  pch = 19, cex = point.size)#add estimates

abline(v=x.axis, lty = 2, lwd = .5, col = "light gray")
axis(1, at = x.tick.marks, label = NA, tick = T)#add tick marks at interval breaks
axis(1, at = x.axis + .5, label = x.label, tick = F, mgp = c(2,.7,0),cex.axis = axis.text)#add labels in middle of intervals
axis(2, at = seq(0,1, by = .1), label = seq(0,100, by = 10), las =2,
    cex.axis = axis.text, mgp = c(2,.7,0))
mtext("Percent Liberal", 2, line = 2.5, cex = axis.text)
mtext("Republican Judges", 3, line = .2, cex = 1.4)
#add legend
top.point.height <- .68
middle.point.height <- top.point.height - .05
bottom.point.height <- top.point.height - .1
x.location <- 6.5
points(x.location, top.point.height, type = "p", col = dd.color,  pch = 19, cex = point.size)
text(x.location + .05, top.point.height, "DD",col = dd.color, adj = 0) 
points(x.location + .1,middle.point.height, type = "p", col = rd.color, cex = point.size, pch = 23, bg = "white")
text(x.location + .15, middle.point.height, "RD",col = rd.color, adj = 0) 
points(x.location + .2,bottom.point.height, type = "p", col = rr.color, cex = point.size, pch = 21)
text(x.location + .25, bottom.point.height, "RR",col = rr.color, adj = 0) 

dev.off()

#The rest of the script converts the info presented in the graphs into tables, which appear in Web Appendix B

#set up exports for tables
#1) aggregate data (by year
#years in rows
  #in columns:  % dem appointees, % of panels with dem maj, 
    #RRR, RRD, DDR, DDD, Mixed

rd <- function (x) {#create shortcut for rounding to 2nd decimal place
  round(x,digits = 2)
  }
  
export.matrix <- cbind(
  x.label,
  paste(rd(ddd.all.ci[,1]),rd(ddd.all.rate),rd(ddd.all.ci[,2]), sep = ", "),
  paste(rd(ddr.all.ci[,1]),rd(ddr.all.rate),rd(ddr.all.ci[,2]), sep = ", "),
  paste(rd(rrd.all.ci[,1]),rd(rrd.all.rate),rd(rrd.all.ci[,2]), sep = ", "),
  paste(rd(rrr.all.ci[,1]),rd(rrr.all.rate),rd(rrr.all.ci[,2]), sep = ", "),
  paste(rd(dem.dd.ci[,1]),rd(dem.dd.rate),rd(dem.dd.ci[,2]), sep = ", "),  
  paste(rd(dem.rd.ci[,1]),rd(dem.rd.rate),rd(dem.rd.ci[,2]), sep = ", "),  
  paste(rd(dem.rr.ci[,1]),rd(dem.rr.rate),rd(dem.rr.ci[,2]), sep = ", "), 
  paste(rd(gop.dd.ci[,1]),rd(gop.dd.rate),rd(gop.dd.ci[,2]), sep = ", "),  
  paste(rd(gop.rd.ci[,1]),rd(gop.rd.rate),rd(gop.rd.ci[,2]), sep = ", "),  
  paste(rd(gop.rr.ci[,1]),rd(gop.rr.rate),rd(gop.rr.ci[,2]), sep = ", ")
)

write.csv(export.matrix, "votes_yearly.table.csv")