#Endemism in avifaunas of Northern Melanesian islands
#Mayr and Diamond (2001, The birds of Northern Melanesia, Oxford University Press) 

library(rgl)


Area_km_2 <- c(49658, 35742, 52, 816, 14, 61, 57, 329, 46, 14, 
		   7174, 1186, 110, 218, 205, 98, 110, 1.5, 414, 41,
		   0.26, 1834, 88, 13, 10, 0.52, 15, 36000, 37, 22512,
		   611, 8591, 8591, 232, 71, 2966, 4095, 381, 5281, 5281,
		   73, 31, 4924, 3566, 109, 544, 2044, 95, 69, 704,
		   704, 805, 640, 142, 13, 35, 502, 381, 122, 4.0,
		   176, 4307, 3090, 42, 11, 20, 684, 20, 9.6, 1.3)

length(Area_km_2)


Distance_km <- c(46, 30, 9, 25, 8, 5, 63, 50, 22, 10,
		     30, 30, 14, 26, 47, 59, 55, 29, 96, 17,
		     74, 241, 43, 235, 170, 174, 162, 171, 63, 51,
		     1, 52, 52, 8, 11, 52, 36, 24, 24, 24,
		     29, 13, 53, 9, 8, 0.6, 0.6, 0.1, 2, 2, 
		     2, 24, 9, 9, 7, 12, 9, 4, 4, 59, 
		     46, 47, 47, 8, 20, 8, 168, 23, 237, 175)

length(Distance_km)
summary(Distance_km)

Endemism_index <- c(1.28, 0.54, 0, 0.04, 0, 0, 0, 0, 0.02, 0,
			  0.41, 0.08, 0.21, 0.10, 0.05, 0.03, 0.05, 0, 0.6, 0,
			  0, 0.72, 0.07, 0.06, 0, 0, 0, 1.81, 0.14, 0.83,
			  0, 0.21, 0.04, 0, 0, 0.10, 0.04, 0.05, 0.29, 0.16,
			  0.02, 0, 0.52, 0.18, 0, 0, 0, 0, 0, 0.13,
			  0, 0.14, 0.09, 0.11, 0, 0.05, 0.08, 0.03, 0.04, 0,
			  0.11, 0.41, 0.88, 0.08, 0.03, 0.02, 0.95, 0, 0, 0)

length(Endemism_index)



plot(log(Area_km_2), Endemism_index)

points(log(Area_km_2[Distance_km > 48.2]), Endemism_index[Distance_km > 48.2], col="red", pch=19)




plot(log(Area_km_2), log(Endemism_index))
points(log(Area_km_2[Distance_km > 48.2]), log(Endemism_index[Distance_km > 48.2]), col="red")


plot(log(Area_km_2), log(Endemism_index), cex=1.5, pch=19, col="gray50")
points(log(Area_km_2[Distance_km > median(Distance_km)]), log(Endemism_index[Distance_km > median(Distance_km)]), col="red", cex=1.5, pch=19)
points(log(Area_km_2[Distance_km <= median(Distance_km)]), log(Endemism_index[Distance_km <= median(Distance_km)]), col="blue", cex=1.5, pch=19)

plot(log(Area_km_2), Endemism_index, cex=1.5, pch=19, col="gray50")
points(log(Area_km_2[Distance_km > median(Distance_km)]), Endemism_index[Distance_km > median(Distance_km)], col="red", cex=1.5, pch=19)
points(log(Area_km_2[Distance_km <= median(Distance_km)]), Endemism_index[Distance_km <= median(Distance_km)], col="blue", cex=1.5, pch=19)

Area_km_2_trim <- Area_km_2[c(1,2,11,13:30,41,43,44,52,56,57,60:65,67:70)]
Distance_km_trim <- Distance_km[c(1,2,11,13:30,41,43,44,52,56,57,60:65,67:70)]
Endemism_index_trim <- Endemism_index[c(1,2,11,13:30,41,43,44,52,56,57,60:65,67:70)]

#par(mar = c(5, 4, 4, 2) + 0.1)
par(mar = c(5, 5, 4, 2) + 0.1)
plot(log(Area_km_2_trim), log(Endemism_index_trim), cex=1.5, pch=19, col="gray40" ,bty="n", xlab=expression(paste("Log (area ", km^2, ")")), xlim=c(2,12), ylim=c(-4,1.5), ylab="Log (endemism index)",  cex.axis=1.5, cex.lab=1.5)
points(log(Area_km_2_trim[Distance_km_trim > median(Distance_km_trim)]), log(Endemism_index_trim[Distance_km_trim > median(Distance_km_trim)]), col="red", cex=1.5, pch=19)
points(log(Area_km_2_trim[Distance_km_trim <= median(Distance_km_trim)]), log(Endemism_index_trim[Distance_km_trim <= median(Distance_km_trim)]), col="blue", cex=1.5, pch=19)
legend(1.8, 1.5, c("< 50 km from colonization source", "> 50 km from colonization source"), pch=19, col=c("red", "blue"), cex=1.5)