################# # looking to see if there is a difference in the ratio of completed to invited reviews from year to year # Molecular Ecology Analysis ---------------------------------------------- # load data ripr.dat <- read.csv("ripr for R.csv") head(ripr.dat) ########### # descriptives ripr.dat$year <- factor(ripr.dat$year) ripr.dat$comp.prop <- ripr.dat$Completed/ripr.dat$Invited ######## #### By Year # proportions tapply(ripr.dat$Completed, ripr.dat$year, sum)/tapply(ripr.dat$Invited, ripr.dat$year, sum) # 2009 2011 2013 2015 # 0.4735915 0.4448052 0.4124629 0.4365325 # use prop.test to get the confidence intervals? tapply(ripr.dat$Completed, ripr.dat$year, sum) # 2009 2011 2013 2015 # 269 274 278 282 tapply(ripr.dat$Invited, ripr.dat$year, sum) # 2009 2011 2013 2015 # 568 616 674 646 # figure: 3 panel: A = invited, B = completed, C = proportion quartz(width = 8, height = 9) nf <- layout(matrix(c(1,2,3,3), 2, 2, byrow = TRUE)) layout.show(nf) par(mar = c(4, 4, 2, 2), cex = 1.1) boxplot(ripr.dat$Invited ~ ripr.dat$year, las = 1, ylab = "Number of Review Invitations", xlab = "Year") text(0.51, 19, "A") par(mar = c(4, 4, 2, 2), cex = 1.1) barplot(table(ripr.dat$year, ripr.dat$Completed), beside = TRUE, las = 1, ylab = "Count", xlab = "Number of Completed Reviews", ylim = c(0, 75), col = c("grey10", "grey45", "grey70", "grey100")); box (bty = "o") legend(18, 72, legend = c("2009", "2011", "2013", "2015"), cex = 0.8, fill = c("grey10", "grey45", "grey70", "grey100"), border = "black") text(1.1, 72, "B") par(mar = c(4, 4, 2, 4), cex = 1.1) plot(c(2009, 2011, 2013, 2015), tapply(ripr.dat$Completed, ripr.dat$year, sum)/tapply(ripr.dat$Invited, ripr.dat$year, sum), ylim = c(0.35, 0.55), las = 1, xlab = "Year", ylab = "proportion of completed reviews", axes = FALSE, pch = 19, type = "b") axis(1, at = c(2009, 2011, 2013, 2015), labels = c("2009", "2011", "2013", "2015")) axis(2, las = 1) box(bty = "o") arrows(2009, 0.4319745, 2009, 0.5155710, angle = 90, code = 3, length = 0) arrows(2011, 0.4052247, 2011, 0.4850839, angle = 90, code = 3, length = 0) arrows(2013, 0.3751689, 2013, 0.4507688, angle = 90, code = 3, length = 0) arrows(2015, 0.4051242, 2015, 0.4837507, angle = 90, code = 3, length = 0) text(2008.85, 0.55, "C") par(new = TRUE) plot(c(2009, 2010, 2011, 2012, 2013, 2014, 2015), c(5.96, 6.457, 5.52, 6.275, 5.84, 6.494, 5.95), axes = FALSE, pch = 3, type = "b", ylab = NA, xlab = NA, ylim = c(5.5, 6.5), lty = 3, col = "grey30") axis(side = 4, las = 1) mtext(side = 4, line = 2.5, "Impact Factor") legend("bottomright", legend = c("proportion completed reviews", "impact factor"), lty = c(1, 3), pch = c(19, 3), col = c("black", "grey30"), cex = 0.8, bty = "n") xx <- barplot(tapply(ripr.dat$Completed, ripr.dat$year, sum)/tapply(ripr.dat$Invited, ripr.dat$year, sum), col = "grey80", ylim = c(0, 0.55), las = 1, ylab = "Mean proportion of completed reviews", width = c(0.5, 0.5, 0.5, 0.5));box(bty = "o") arrows(xx[1], 0.4319745, xx[1], 0.5155710, angle = 90, code = 3) arrows(xx[2], 0.4052247, xx[2], 0.4850839, angle = 90, code = 3) arrows(xx[3], 0.3751689, xx[3], 0.4507688, angle = 90, code = 3) arrows(xx[4], 0.4051242, xx[4], 0.4837507, angle = 90, code = 3) text(0.05, 0.53, "C") # calculating confidence intervals for the figure above prop.test(269, 568) # 0.4319745 0.5155710 prop.test(274, 616) # 0.4052247 0.4850839 prop.test(278, 674) # 0.3751689 0.4507688 prop.test(282, 635) # 0.4051242 0.4837507 ###### By Editorial decision category # Decision # r, re, rar, mr levels(ripr.dat$Decision) <- c("MRev", "Reject", "RAR", "REnc") levels(ripr.dat$Decision) ripr.dat$Decision <- relevel(ripr.dat$Decision, ref = "MRev") ripr.dat$Decision <- relevel(ripr.dat$Decision, ref = "RAR") ripr.dat$Decision <- relevel(ripr.dat$Decision, ref = "REnc") ripr.dat$Decision <- relevel(ripr.dat$Decision, ref = "Reject") # make figure quartz(width = 8, height = 9) nf <- layout(matrix(c(1,2,3,3), 2, 2, byrow = TRUE)) layout.show(nf) par(mar = c(4, 4, 2, 2), cex = 1.1) boxplot(ripr.dat$Invited ~ ripr.dat$Decision, las = 1, ylab = "Number of Review Invitations", names = FALSE) text(0.51, 19, "A") mtext(levels(ripr.dat$Decision), at = c(1,2,3,4), side = 1, line = 1, cex = 0.9) mtext("Decision", side = 1, line = 2) par(mar = c(4, 4, 2, 2), cex = 1.1) barplot(table(ripr.dat$Decision, ripr.dat$Completed), beside = TRUE, las = 1, ylab = "Count", xlab = "Number of Completed Reviews", ylim = c(0, 130), col = c("grey10", "grey45", "grey70", "grey100")); box (bty = "o") legend(16, 125, legend = levels(ripr.dat$Decision), cex = 0.8, fill = c("grey10", "grey45", "grey70", "grey100"), border = "black") text(1.1, 125, "B") par(mar = c(4, 4, 2, 2), cex = 1.1) xx <- c(1,2,3,4) plot(xx, tapply(ripr.dat$Completed, ripr.dat$Decision, sum)/tapply(ripr.dat$Invited, ripr.dat$Decision, sum), ylim = c(0.33, 0.55), las = 1, xlab = "Decision", ylab = "proportion of completed reviews", axes = FALSE, pch = 19, type = "b") axis(1, at = xx, labels = c("Reject", "REnc", "RAR", "MRev")) axis(2, las = 1) box(bty = "o") arrows(xx[1], 0.3853596, xx[1], 0.4486490, angle = 90, code = 3, length = 0) arrows(xx[2], 0.4382238, xx[2], 0.5003070, angle = 90, code = 3, length = 0) arrows(xx[3], 0.3923736, xx[3], 0.4964988, angle = 90, code = 3, length = 0) arrows(xx[4], 0.3416978, xx[4], 0.5059918, angle = 90, code = 3, length = 0) text(0.93, 0.55, "C") xx <- barplot(tapply(ripr.dat$Completed, ripr.dat$Decision, sum)/tapply(ripr.dat$Invited, ripr.dat$Decision, sum), col = "grey80", ylim = c(0, 0.55), las = 1, ylab = "Mean proportion of completed reviews", width = c(0.5, 0.5, 0.5, 0.5));box(bty = "o") arrows(xx[1], 0.3853596, xx[1], 0.4486490, angle = 90, code = 3) arrows(xx[2], 0.4382238, xx[2], 0.5003070, angle = 90, code = 3) arrows(xx[3], 0.3923736, xx[3], 0.4964988, angle = 90, code = 3) arrows(xx[4], 0.3416978, xx[4], 0.5059918, angle = 90, code = 3) text(0.05, 0.53, "C") ##### # calculate cinfidence intervals for the figure above # so the cis are: prop.test(400, 960) # 0.3853596 0.4486490 prop.test(479, 1021) # 0.4382238 0.5003070 prop.test(162, 365) # 0.3923736 0.4964988 prop.test(62, 147) # 0.3416978 0.5059918 ########## # GLM of proportion completed by year m5 <- glm(cbind(Completed, Invited..includes.agree.decline.completed.Invalid.Email..excludes.Uninvited. - Completed) ~ year, family = "binomial", data = ripr.dat) summary(m5) # Call: # glm(formula = cbind(Completed, Invited..includes.agree.decline.completed.Invalid.Email..excludes.Uninvited. - # Completed) ~ year, family = "binomial", data = ripr.dat) # # Deviance Residuals: # Min 1Q Median 3Q Max # -2.6799 -0.5566 0.2215 0.8871 2.6617 # # Coefficients: # Estimate Std. Error z value Pr(>|z|) # (Intercept) -0.10573 0.08404 -1.258 0.2083 # year2011 -0.11595 0.11677 -0.993 0.3207 # year2013 -0.24806 0.11482 -2.160 0.0307 * # year2015 -0.14951 0.11557 -1.294 0.1957 # --- # Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 # # (Dispersion parameter for binomial family taken to be 1) # # Null deviance: 469.02 on 399 degrees of freedom # Residual deviance: 464.26 on 396 degrees of freedom # AIC: 1300.3 ## LRT for p-value drop1(m5, test = "Chis") ## OR and CI exp(cbind(m5$coef, confint(m5))) ########## # GLM of proportion completed by decision m6 <- glm(cbind(Completed, Invited..includes.agree.decline.completed.Invalid.Email..excludes.Uninvited. - Completed) ~ year + Decision, family = "binomial", data = ripr.dat) summary(m6) Call: # glm(formula = cbind(Completed, Invited..includes.agree.decline.completed.Invalid.Email..excludes.Uninvited. - # Completed) ~ year + Decision, family = "binomial", data = ripr.dat) # # Deviance Residuals: # Min 1Q Median 3Q Max # -2.7272 -0.5272 0.2033 0.9820 2.6825 # # Coefficients: # Estimate Std. Error z value Pr(>|z|) # (Intercept) -0.19000 0.18302 -1.038 0.299 # year2011 -0.11042 0.11718 -0.942 0.346 # year2013 -0.26836 0.11533 -2.327 0.020 * # year2015 -0.16972 0.11601 -1.463 0.143 # Decisionr -0.01759 0.17983 -0.098 0.922 # Decisionrar 0.09699 0.19780 0.490 0.624 # Decisionre 0.20861 0.17944 1.163 0.245 # --- # Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 # # (Dispersion parameter for binomial family taken to be 1) # # Null deviance: 469.02 on 399 degrees of freedom # Residual deviance: 457.77 on 393 degrees of freedom # AIC: 1299.8 #LRT for p-value drop1(m6, test = "Chis") # OR and CI exp(cbind(m6$coef, confint(m6))) ########## #means for Table 1 tapply(ripr.dat$Declined, ripr.dat$year, sum)/tapply(ripr.dat$Declined, ripr.dat$year, length) # 2009 2011 2013 2015 # 2.76 3.18 3.77 3.43 tapply(ripr.dat$Unassigned, ripr.dat$year, sum)/tapply(ripr.dat$Unassigned, ripr.dat$year, length) # 2009 2011 2013 2015 # 0.23 0.24 0.19 0.10 tapply(ripr.dat$Uninvited, ripr.dat$year, sum)/tapply(ripr.dat$Uninvited, ripr.dat$year, length) # 2009 2011 2013 2015 # 0.14 0.08 0.05 0.02 # proportions for decision tapply(ripr.dat$Completed, ripr.dat$Decision, sum)/tapply(ripr.dat$Invited, ripr.dat$Decision, sum) # mr r rar re # 0.4217687 0.4153686 0.4402174 0.4668616 # use prop.test to get the confidence intervals? tapply(ripr.dat$Completed, ripr.dat$Decision, sum) # Reject REnc RAR MRev # 400 479 162 62 tapply(ripr.dat$Completed, ripr.dat$Decision, sum)/tapply(ripr.dat$Completed, ripr.dat$Decision, length) # mr r rar re # 2.695652 2.666667 2.745763 2.851190 tapply(ripr.dat$Invited, ripr.dat$Decision, sum) # mr r rar re # 147 963 368 1026 tapply(ripr.dat$Invited, ripr.dat$Decision, sum)/tapply(ripr.dat$Invited, ripr.dat$Decision, length) # mr r rar re # 6.391304 6.420000 6.237288 6.107143 # BES data analysis ------------------------------------------------------- ######################## # additional analysis with BES data ######################### load the data jpe.dat <- read.csv("~/Google Drive/RIPR paper/Petchey et al 2014 PloS ONE data/alldatafiles/authors/JPE authors.csv") head(jpe.dat) jec.dat <- read.csv("~/Google Drive/RIPR paper/Petchey et al 2014 PloS ONE data/alldatafiles/authors/JEC authors.csv") head(jec.dat) jae.dat <- read.csv("~/Google Drive/RIPR paper/Petchey et al 2014 PloS ONE data/alldatafiles/authors/JAE authors.csv") head(jae.dat) fec.dat <- read.csv("~/Google Drive/RIPR paper/Petchey et al 2014 PloS ONE data/alldatafiles/authors/FEC authors.csv") head(fec.dat) # jpe.dat <- read.csv("JPE authors.csv") # head(jpe.dat) # jec.dat <- read.csv("JEC authors.csv") # head(jec.dat) # jae.dat <- read.csv("JAE authors.csv") # head(jae.dat) # fec.dat <- read.csv("FEC authors.csv") # head(fec.dat) ########## # JPE ########## ### analyze journals separately # remove papers with n=0 invited jpe.dat <- jpe.dat[-which(jpe.dat$X..Reviewers.Invited == 0), ] levels(factor(jpe.dat$Submission.Year)) jpe.dat[which(jpe.dat$X..Reviews.Completed == 0), ] jpe.dat$comp.prop <- jpe.dat$X..Reviews.Completed/jpe.dat$X..Reviewers.Invited hist(jpe.dat$comp.prop) # one is greater than 1, so remove jpe.dat <- jpe.dat[-which(jpe.dat$comp.prop >1), ] # what do the decisions look like? levels(jpe.dat$Decision) jpe.dat$Decision <- droplevels(jpe.dat$Decision) levels(jpe.dat$Decision) # remove the ones with immediate rejections: jpe.dat[which(jpe.dat$Decision == "Immediate Reject - Executive Editor"), ] jpe.dat <- jpe.dat[-which(jpe.dat$Decision == "Immediate Reject - Executive Editor"), ] ## GLM by year: j1 <- glm(cbind(X..Reviews.Completed, X..Reviewers.Invited - X..Reviews.Completed) ~ Submission.Year, data = jpe.dat, family = "binomial") summary(j1) # Call: # glm(formula = cbind(X..Reviews.Completed, X..Reviewers.Invited - # X..Reviews.Completed) ~ Submission.Year, family = "binomial", # data = jpe.dat) # # Deviance Residuals: # Min 1Q Median 3Q Max # -3.9496 -0.5034 0.1958 0.9154 2.5234 # # Coefficients: # Estimate Std. Error z value Pr(>|z|) # (Intercept) 122.932077 14.603003 8.418 <2e-16 *** # Submission.Year -0.061258 0.007276 -8.419 <2e-16 *** # --- # Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 # # (Dispersion parameter for binomial family taken to be 1) # # Null deviance: 4326.4 on 3517 degrees of freedom # Residual deviance: 4255.2 on 3516 degrees of freedom # AIC: 9947.4 drop1(j1, test = "Chi") # cbind(X..Reviews.Completed, X..Reviewers.Invited - X..Reviews.Completed) ~ # Submission.Year # Df Deviance AIC LRT Pr(>Chi) # 4255.2 9947.4 # Submission.Year 1 4326.4 10016.6 71.17 < 2.2e-16 *** exp(cbind(j1$coef, confint(j1))) # 2.5 % 97.5 % # (Intercept) 2.447501e+53 9.249933e+40 6.719268e+65 # Submission.Year 9.405808e-01 9.272530e-01 9.540827e-01 # The coefficient suggests that for every increase in 1 year in the dataset, the odds of a completed review goes down by 6% # predicted probabilities newdata <- data.frame(Submission.Year = 2003) predict(j1, newdata, type = "response") mean(jpe.dat$comp.prop[jpe.dat$Submission.Year == 2003]) newdata <- data.frame(Submission.Year = 2010) predict(j1, newdata, type = "response") mean(jpe.dat$comp.prop[jpe.dat$Submission.Year == 2010]) # get the averages by year tapply(jpe.dat$X..Reviews.Completed, jpe.dat$Submission.Year, sum) # 2003 2004 2005 2006 2007 2008 2009 2010 # 618 862 1068 1190 1138 969 1075 1132 tapply(jpe.dat$X..Reviewers.Invited, jpe.dat$Submission.Year, sum) # 2003 2004 2005 2006 2007 2008 2009 2010 # 1129 1606 1965 2285 2338 1995 2385 2456 # how many papers going out for review each year? summary(factor(jpe.dat$Submission.Year)) # 2003 2004 2005 2006 2007 2008 2009 2010 # 255 363 446 507 502 420 488 537 tapply(jpe.dat$X..Reviews.Completed, jpe.dat$Submission.Year, sum)/tapply(jpe.dat$X..Reviewers.Invited, jpe.dat$Submission.Year, sum) # 2003 2004 2005 2006 2007 2008 2009 2010 # 0.5473871 0.5367372 0.5435115 0.5207877 0.4867408 0.4857143 0.4507338 0.4609121 # confidence intervals prop.test(618, 1129) # 0.5177971 0.5766508 prop.test(862, 1606) # 0.0.5119792 0.561317 prop.test(1068, 1965) # 0.5211695 0.5656817 prop.test(1190, 2285) # 0.5000679 0.5414370 prop.test(1138, 2338) # 0.4663057 0.5072198 prop.test(969, 1995) # 0.4635814 0.5079027 prop.test(1075, 2385) # 0.4306515 0.4709762 prop.test(1132, 2456) # 0.4410719 0.4808756 # plot the results xv.jpe <- seq(2003, 2010, 0.1) plot(jpe.dat$Submission.Year, jpe.dat$comp.prop, ylab = "Proportion Complete", jitter = "TRUE") lines(xv.jpe, predict(j1, list(Submission.Year = xv.jpe), type = "response")) ################ # JEC ################ # remove n=0 invitations: jec.dat <- jec.dat[-which(jec.dat $X..Reviewers.Invited == 0), ] levels(factor(jec.dat$Submission.Year)) jec.dat$comp.prop <- jec.dat$X..Reviews.Completed/jec.dat$X..Reviewers.Invited hist(jec.dat$comp.prop) # chicking decision categories levels(jec.dat$Decision) jec.dat $Decision <- droplevels(jec.dat $Decision) levels(jec.dat$Decision) table(jec.dat$Decision, jec.dat$Submission.Year) #GLM with year j2 <- glm(cbind(X..Reviews.Completed, X..Reviewers.Invited - X..Reviews.Completed) ~ Submission.Year, data = jec.dat, family = "binomial") summary(j2) # glm(formula = cbind(X..Reviews.Completed, X..Reviewers.Invited - # X..Reviews.Completed) ~ Submission.Year, family = "binomial", # data = jec.dat) # # Deviance Residuals: # Min 1Q Median 3Q Max # -3.2017 -0.4862 0.2927 1.4084 2.2440 # # Coefficients: # Estimate Std. Error z value Pr(>|z|) # (Intercept) 205.542398 16.240922 12.66 <2e-16 *** # Submission.Year -0.102396 0.008092 -12.65 <2e-16 *** # --- # Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 # # (Dispersion parameter for binomial family taken to be 1) # # Null deviance: 3890.2 on 3093 degrees of freedom # Residual deviance: 3728.1 on 3092 degrees of freedom # AIC: 8050.4 drop1(j2, test = "Chi") # cbind(X..Reviews.Completed, X..Reviewers.Invited - X..Reviews.Completed) ~ # Submission.Year # Df Deviance AIC LRT Pr(>Chi) # 3728.1 8050.4 # Submission.Year 1 3890.2 8210.5 162.12 < 2.2e-16 *** exp(cbind(j2$coef, confint(j2))) # 2.5 % 97.5 % # (Intercept) 1.844715e+89 2.873284e+75 1.283811e+103 # Submission.Year 9.026720e-01 8.884502e-01 9.170846e-01 # for every year, 10% reduction in the odds of a completed review # get the averages by year tapply(jec.dat$X..Reviews.Completed, jec.dat$Submission.Year, sum) # 2003 2004 2005 2006 2007 2008 2009 2010 # 553 727 787 769 794 711 842 1041 tapply(jec.dat$X..Reviewers.Invited, jec.dat$Submission.Year, sum) # 2003 2004 2005 2006 2007 2008 2009 2010 # 908 1198 1370 1577 1560 1480 1773 2417 tapply(jec.dat$X..Reviews.Completed, jec.dat$Submission.Year, sum)/tapply(jec.dat$X..Reviewers.Invited, jec.dat$Submission.Year, sum) # 2003 2004 2005 2006 2007 2008 2009 2010 # 0.6090308 0.6068447 0.5744526 0.4876347 0.5089744 0.4804054 0.4749013 0.4306992 prop.test(553, 908) # 0.5763396 0.6407875 prop.test(727, 1198) # 0.5784649 0.6345312 prop.test(787, 1370) # 0.5477322 0.6007508 prop.test(769, 1577) # 0.4627089 0.5126215 prop.test(794, 1560) # 0.4838551 0.5340489 prop.test(711, 1480) # 0.4546989 0.5062147 prop.test(842, 1773) # 0.4514556 0.4984568 prop.test(1041, 2417) # 0.4108782 0.4507424 # plot the results xv.jec ########### #JAE ########### # remove n = 0 invites jae.dat <- jae.dat[-which(jae.dat $X..Reviewers.Invited == 0), ] levels(factor(jae.dat $Submission.Year)) jae.dat $comp.prop <- jae.dat $X..Reviews.Completed/jae.dat $X..Reviewers.Invited hist(jae.dat $comp.prop) levels(jae.dat $Decision) jae.dat $Decision <- droplevels(jae.dat $Decision) levels(jae.dat $Decision) # remove immediate reject: jae.dat[which(jae.dat$Decision == "Return - without review"), ] table(jae.dat$X..Reviews.Completed, jae.dat$Submission.Year) ##GLM against year j3 <- glm(cbind(X..Reviews.Completed, X..Reviewers.Invited - X..Reviews.Completed) ~ Submission.Year, data = jae.dat, family = "binomial") summary(j3) # Deviance Residuals: # Min 1Q Median 3Q Max # -4.0582 -0.4905 0.3317 1.4629 2.3666 # # Coefficients: # Estimate Std. Error z value Pr(>|z|) # (Intercept) 103.351520 15.723399 6.573 4.93e-11 *** # Submission.Year -0.051426 0.007835 -6.563 5.26e-11 *** # --- # Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 # # (Dispersion parameter for binomial family taken to be 1) # # Null deviance: 4211.2 on 3440 degrees of freedom # Residual deviance: 4168.0 on 3439 degrees of freedom # AIC: 8708 drop1(j3, test = "Chis") # cbind(X..Reviews.Completed, X..Reviewers.Invited - X..Reviews.Completed) ~ # Submission.Year # Df Deviance AIC LRT Pr(>Chi) # 4168.0 8708.0 # Submission.Year 1 4211.2 8749.2 43.224 4.881e-11 *** exp(cbind(j3$coef, confint(j3))) # 2.5 % 97.5 % # (Intercept) 7.673526e+44 3.231485e+31 1.897833e+58 # Submission.Year 9.498743e-01 9.353887e-01 9.645648e-01 # reduction in Odds of 5% per year # get the averages by year tapply(jae.dat$X..Reviews.Completed, jae.dat$Submission.Year, sum) # 2003 2004 2005 2006 2007 2008 2009 2010 # 698 855 761 957 987 869 847 890 tapply(jae.dat$X..Reviewers.Invited, jae.dat$Submission.Year, sum) # 2003 2004 2005 2006 2007 2008 2009 2010 # 1257 1547 1361 1679 1663 1593 1758 1896 tapply(jae.dat$X..Reviews.Completed, jae.dat$Submission.Year, sum)/tapply(jae.dat$X..Reviewers.Invited, jae.dat$Submission.Year, sum) # 2003 2004 2005 2006 2007 2008 2009 2010 # 0.5552904 0.5526826 0.5591477 0.5699821 0.5935057 0.5455116 0.4817975 0.4694093 prop.test(698, 1257) # 0.5272930 0.5829459 prop.test(855, 1547) # 0.5274814 0.5776195 prop.test(761, 1361) # 0.5322720 0.5856858 prop.test(957, 1679) # 0.5458690 0.5937718 prop.test(987, 1663) # 0.5694067 0.6171683 prop.test(869, 1593) # 0.5206654 0.5701360 prop.test(847, 1758) # 0.4582221 0.5054532 prop.test(890, 1896) # 0.4467675 0.4921763 tapply(predict(j3, type = "response"), jae.dat$Submission.Year, mean) ######## #FEC ######## # remove n = 0 invites fec.dat <- fec.dat[-which(fec.dat $X..Reviewers.Invited == 0), ] levels(factor(fec.dat $Submission.Year)) fec.dat $comp.prop <- fec.dat $X..Reviews.Completed/fec.dat $X..Reviewers.Invited hist(fec.dat $comp.prop) levels(fec.dat $Decision) fec.dat $Decision <- droplevels(fec.dat $Decision) levels(fec.dat $Decision) # remove immediate reject: fec.dat[which(fec.dat $Decision == "Not appropriate - immediate"), ] ## GLM against year j4 <- glm(cbind(X..Reviews.Completed, X..Reviewers.Invited - X..Reviews.Completed) ~ Submission.Year, data = fec.dat, family = "binomial") summary(j4) # Deviance Residuals: # Min 1Q Median 3Q Max # -3.7153 -0.5063 0.4388 1.1250 2.4439 # # Coefficients: # Estimate Std. Error z value Pr(>|z|) # (Intercept) 78.21052 15.63319 5.003 5.65e-07 *** # Submission.Year -0.03896 0.00779 -5.002 5.69e-07 *** # --- # Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 # # (Dispersion parameter for binomial family taken to be 1) # # Null deviance: 3730.8 on 3012 degrees of freedom # Residual deviance: 3705.8 on 3011 degrees of freedom # AIC: 8090.8 exp(cbind(j4$coef, confint(j4))) # 2.5 % 97.5 % # (Intercept) 9.255463e+33 4.625514e+20 1.907191e+47 # Submission.Year 9.617867e-01 9.472059e-01 9.765778e-01 # odds reduced by 4% per year # get the averages by year tapply(fec.dat$X..Reviews.Completed, fec.dat$Submission.Year, sum) # 2003 2004 2005 2006 2007 2008 2009 2010 # 669 723 716 777 807 797 780 954 tapply(fec.dat$X..Reviewers.Invited, fec.dat$Submission.Year, sum) # 2003 2004 2005 2006 2007 2008 2009 2010 # 1309 1312 1414 1441 1561 1535 1703 2047 tapply(fec.dat$X..Reviews.Completed, fec.dat$Submission.Year, sum)/tapply(fec.dat$X..Reviewers.Invited, fec.dat$Submission.Year, sum) # 2003 2004 2005 2006 2007 2008 2009 2010 # 0.5110772 0.5510671 0.5063649 0.5392089 0.5169763 0.5192182 0.4580153 0.4660479 ########### # Fox et al data for FEC ########### fec.fox.dat <- read.csv("~/Google Drive/RIPR paper/Fox et al 2016 Functional Ecology data/Fec for R.csv") tapply(fec.fox.dat$agreed, fec.fox.dat$year, sum)/tapply(fec.fox.dat$invited, fec.fox.dat $year, sum) # 2004 2005 2006 2007 2008 2009 2010 2011 2012 2013 2014 # 0.6390698 0.6682199 0.6518456 0.5447552 0.5386799 0.4777487 0.4807796 0.4224757 0.4104355 0.4053816 0.3853142 # Figure S1 quartz(width = 8, height = 5) par(cex=1.1, mar = c(5, 4, 3, 2)) plot(seq(2003, 2010), tapply(fec.dat$X..Reviews.Completed, fec.dat$Submission.Year, sum)/tapply(fec.dat$X..Reviewers.Invited, fec.dat$Submission.Year, sum), xlim = c(2003, 2014), pch = 16, type = "b", las = 1, ylim = c(0.35, 0.70), ylab = "Proportion completed or agreed reviews", xlab = "Year") # arrows(2003, 0.4836235, 2003, 0.5384651, angle = 90, code = 3, length = 0) # arrows(2004, 0.5236611, 2004, 0.5781707, angle = 90, code = 3, length = 0) # arrows(2005, 0.4799711, 2005, 0.5327238, angle = 90, code = 3, length = 0) # arrows(2006, 0.5130549, 2006, 0.5651516, angle = 90, code = 3, length = 0) # arrows(2007, 0.4918556, 2007, 0.5420126, angle = 90, code = 3, length = 0) # arrows(2008, 0.4938815, 2008, 0.5444578, angle = 90, code = 3, length = 0) # arrows(2009, 0.4341812, 2009, 0.4820407, angle = 90, code = 3, length = 0) # arrows(2010, 0.4442785, 2010, 0.4879459, angle = 90, code = 3, length = 0) points(seq(2004, 2014), tapply(fec.fox.dat$agreed, fec.fox.dat$year, sum)/tapply(fec.fox.dat$invited, fec.fox.dat $year, sum), pch = 1, type = "b") legend("topright", pch = c(1, 16), legend = c("Fox et al 2016a", "Petchey et al 2014a")) # calculate CIs for figure prop.test(tapply(fec.fox.dat$agreed, fec.fox.dat$year, sum)[1], tapply(fec.fox.dat$invited, fec.fox.dat$year, sum)[1]) # 0.6094402 0.6676929 prop.test(tapply(fec.fox.dat$agreed, fec.fox.dat$year, sum)[2], tapply(fec.fox.dat$invited, fec.fox.dat$year, sum)[2]) # 0.6390165 0.6962034 prop.test(tapply(fec.fox.dat$agreed, fec.fox.dat$year, sum)[3], tapply(fec.fox.dat$invited, fec.fox.dat$year, sum)[3]) # 0.6239274 0.6787732 prop.test(tapply(fec.fox.dat$agreed, fec.fox.dat$year, sum)[4], tapply(fec.fox.dat$invited, fec.fox.dat$year, sum)[4]) # 0.5185084 0.5707590 prop.test(tapply(fec.fox.dat$agreed, fec.fox.dat$year, sum)[5], tapply(fec.fox.dat$invited, fec.fox.dat$year, sum)[5]) # 0.5122255 0.5649211 prop.test(tapply(fec.fox.dat$agreed, fec.fox.dat$year, sum)[6], tapply(fec.fox.dat$invited, fec.fox.dat$year, sum)[6]) # 0.4524649 0.5031455 prop.test(tapply(fec.fox.dat$agreed, fec.fox.dat$year, sum)[7], tapply(fec.fox.dat$invited, fec.fox.dat$year, sum)[7]) # 0.4577877 0.5038523 prop.test(tapply(fec.fox.dat$agreed, fec.fox.dat$year, sum)[8], tapply(fec.fox.dat$invited, fec.fox.dat$year, sum)[8]) # 0.4025083 0.4426970 prop.test(tapply(fec.fox.dat$agreed, fec.fox.dat$year, sum)[9], tapply(fec.fox.dat$invited, fec.fox.dat$year, sum)[9]) # 0.3908456 0.4303106 prop.test(tapply(fec.fox.dat$agreed, fec.fox.dat$year, sum)[10], tapply(fec.fox.dat$invited, fec.fox.dat$year, sum)[10]) # 0.3851293 0.4259574 prop.test(tapply(fec.fox.dat$agreed, fec.fox.dat$year, sum)[11], tapply(fec.fox.dat$invited, fec.fox.dat$year, sum)[11]) # 0.359066 0.412237 prop.test(669, 1309) # 0.4836235 0.5384651 prop.test(723, 1312) # 0.5236611 0.5781707 prop.test(716, 1414) # 0.4799711 0.5327238 prop.test(777, 1441) # 0.5130549 0.5651516 prop.test(807, 1561) # 0.4918556 0.5420126 prop.test(797, 1535) # 0.4938815 0.5444578 prop.test(780, 1703) # 0.4341812 0.4820407 prop.test(954, 2047) # 0.4442785 0.4879459 ############## # figures for the Petchey et al BES data #AUGUST 19TH: ADD IMPACT FACTORS TO THE RIGHT HAND AXIS jpe.if <- c(3.205, 3.266, 4.594, 4.527, 4.220, 4.560, 4.197, 4.970) jec.if <- c(2.833, 3.397, 4.277, 4.239, 4.422, 4.262, 4.690, 5.260) jae.if <- c(2.843, 3.342, 3.399, 3.390, 3.747, 4.220, 3.714, 4.457) fec.if <- c(2.351, 2.718, 3.149, 3.417, 3.157, 3.699, 4.546, 4.645) quartz(width = 8, height = 12) nf <- layout(matrix(c(1,2,3,4), 4, 1, byrow = TRUE)) layout.show(nf) par(mar = c(2, 4, 1, 4), cex = 1.1) plot(seq(1, 8), as.vector(tapply(jpe.dat$X..Reviews.Completed, jpe.dat$Submission.Year, sum)/tapply(jpe.dat$X..Reviewers.Invited, jpe.dat$Submission.Year, sum)), pch = 19, ylim = c(0.4, 0.7), las = 1, ylab = "Proportion completed reviews", xlab = "Year", axes = FALSE, type = "b");box(bty = "o") axis(1, at = seq(1, 8), labels = seq(2003, 2010)) axis(2, las = 1) arrows(1, 0.52, 1, 0.58, angle = 90, code = 3, length = 0) arrows(2, 0.51, 2, 0.56, angle = 90, code = 3, length = 0) arrows(3, 0.52, 3, 0.57, angle = 90, code = 3, length = 0) arrows(4, 0.50, 4, 0.54, angle = 90, code = 3, length = 0) arrows(5, 0.47, 5, 0.51, angle = 90, code = 3, length = 0) arrows(6, 0.46, 6, 0.51, angle = 90, code = 3, length = 0) arrows(7, 0.43, 7, 0.47, angle = 90, code = 3, length = 0) arrows(8, 0.44, 8, 0.48, angle = 90, code = 3, length = 0) text(1, 0.67, "A: Journal of Applied Ecology", adj = 0) par(new = TRUE) plot(seq(1, 8), jpe.if, axes = FALSE, pch = 3, type = "b", ylab = NA, ylim = c(2.3, 5.3), lty = 3, col = "grey30") axis(side = 4, las = 1) mtext(side = 4, line = 2.5, "Impact Factor") legend(3, 2.8, legend = c("proportion completed reviews", "impact factor"), lty = c(1, 3), pch = c(19, 3), col = c("black", "grey30"), cex = 0.8, bty = "n") par(mar = c(2, 4, 1, 4), cex = 1.1) plot(seq(1, 8), as.vector(tapply(jec.dat$X..Reviews.Completed, jec.dat $Submission.Year, sum)/tapply(jec.dat $X..Reviewers.Invited, jec.dat $Submission.Year, sum)), pch = 19, ylim = c(0.4, 0.7), las = 1, ylab = "Proportion completed reviews", xlab = "Year", axes = FALSE, type = "b");box(bty = "o") axis(1, at = seq(1, 8), labels = seq(2003, 2010)) axis(2, las = 1) arrows(1, 0.5763396, 1, 0.6407875, angle = 90, code = 3, length = 0) arrows(2, 0.5784649, 2, 0.6345312, angle = 90, code = 3, length = 0) arrows(3, 0.5477322, 3, 0.6007508, angle = 90, code = 3, length = 0) arrows(4, 0.4627089, 4, 0.5126215, angle = 90, code = 3, length = 0) arrows(5, 0.4838551, 5, 0.5340489, angle = 90, code = 3, length = 0) arrows(6, 0.4546989, 6, 0.5062147, angle = 90, code = 3, length = 0) arrows(7, 0.4514556, 7, 0.4984568, angle = 90, code = 3, length = 0) arrows(8, 0.4108782, 8, 0.4507424, angle = 90, code = 3, length = 0) text(1, 0.67, "B: Journal of Ecology", adj = 0) par(new = TRUE) plot(seq(1, 8), jec.if, axes = FALSE, pch = 3, type = "b", ylab = NA, ylim = c(2.3, 5.3), lty = 3, col = "grey30") axis(side = 4, las = 1) mtext(side = 4, line = 2.5, "Impact Factor") legend(3, 2.8, legend = c("proportion completed reviews", "impact factor"), lty = c(1, 3), pch = c(19, 3), col = c("black", "grey30"), cex = 0.8, bty = "n") par(mar = c(2, 4, 1, 4), cex = 1.1) plot(seq(1, 8), as.vector(tapply(jae.dat$X..Reviews.Completed, jae.dat $Submission.Year, sum)/tapply(jae.dat $X..Reviewers.Invited, jae.dat $Submission.Year, sum)), pch = 19, ylim = c(0.4, 0.7), las = 1, ylab = "Proportion completed reviews", xlab = "Year", axes = FALSE, type = "b");box(bty = "o") axis(1, at = seq(1, 8), labels = seq(2003, 2010)) axis(2, las = 1) arrows(1, 0.5272930, 1, 0.5829459, angle = 90, code = 3, length = 0) arrows(2, 0.5274814, 2, 0.5776195, angle = 90, code = 3, length = 0) arrows(3, 0.5322720, 3, 0.5856858, angle = 90, code = 3, length = 0) arrows(4, 0.5458690, 4, 0.5937718, angle = 90, code = 3, length = 0) arrows(5, 0.5694067, 5, 0.6171683, angle = 90, code = 3, length = 0) arrows(6, 0.5206654, 6, 0.5701360, angle = 90, code = 3, length = 0) arrows(7, 0.4582221, 7, 0.5054532, angle = 90, code = 3, length = 0) arrows(8, 0.4467675, 8, 0.4921763, angle = 90, code = 3, length = 0) text(1, 0.67, "C: Journal of Animal Ecology", adj = 0) par(new = TRUE) plot(seq(1, 8), jae.if, axes = FALSE, pch = 3, type = "b", ylab = NA, ylim = c(2.3, 5.3), lty = 3, col = "grey30") axis(side = 4, las = 1) mtext(side = 4, line = 2.5, "Impact Factor") legend(3, 2.8, legend = c("proportion completed reviews", "impact factor"), lty = c(1, 3), pch = c(19, 3), col = c("black", "grey30"), cex = 0.8, bty = "n") par(mar = c(2, 4, 1, 4), cex = 1.1) plot(seq(1, 8), as.vector(tapply(fec.dat$X..Reviews.Completed, fec.dat $Submission.Year, sum)/tapply(fec.dat $X..Reviewers.Invited, fec.dat $Submission.Year, sum)), pch = 19, ylim = c(0.4, 0.7), las = 1, ylab = "Proportion completed reviews", xlab = "Year", axes = FALSE, type = "b");box(bty = "o") axis(1, at = seq(1, 8), labels = seq(2003, 2010)) axis(2, las = 1) arrows(1, 0.4836235, 1, 0.5384651, angle = 90, code = 3, length = 0) arrows(2, 0.5236611, 2, 0.5781707, angle = 90, code = 3, length = 0) arrows(3, 0.4799711, 3, 0.5327238, angle = 90, code = 3, length = 0) arrows(4, 0.5130549, 4, 0.5651516, angle = 90, code = 3, length = 0) arrows(5, 0.4918556, 5, 0.5420126, angle = 90, code = 3, length = 0) arrows(6, 0.4938815, 6, 0.5444578, angle = 90, code = 3, length = 0) arrows(7, 0.4341812, 7, 0.4820407, angle = 90, code = 3, length = 0) arrows(8, 0.4442785, 8, 0.4879459, angle = 90, code = 3, length = 0) text(1, 0.67, "D: Functional Ecology", adj = 0) par(new = TRUE) plot(seq(1, 8), fec.if, axes = FALSE, pch = 3, type = "b", ylab = NA, ylim = c(2.3, 5.3), lty = 3, col = "grey30") axis(side = 4, las = 1) mtext(side = 4, line = 2.5, "Impact Factor") legend(3, 2.8, legend = c("proportion completed reviews", "impact factor"), lty = c(1, 3), pch = c(19, 3), col = c("black", "grey30"), cex = 0.8, bty = "n")