Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.
This week we have a guest blog post by Neil Kelley
Last week, I stumbled on Vijay’s blog post demonstrating his new package rvertnet. Although I am a paleontologist, some of my research involves anatomical comparison between extinct species and extant relatives or ecological analogs, so I have some experience using VertNet to track down specimens of interest in museum collections around the country.
Although I have been using R off an on for years, I have always been a little terrified of it. This has less to do with R–which is really quite intuitive and forgiving once you get the hang of it–and more to do with my phobia of comand-line interfaces. I’ve been spoiled by rich GUIs, dropdown menus and toolbars. The ominous pulse of a cursor in an empty console window can still send me into a cold sweat.
But lately that is starting to change thanks to my discovery of a number of resources including RStudio, sites like Quick-R and R Cookbook, the ever growing collection of valuable discussion threads on StackOverflow and the R-help mailing list, and of course the proliferation of R-themed blogs like Vijay’s. I suppose the true sign that I am overcoming my fear is that I find myself “playing” with R to learn and practice skills that I can use in my own research.
So I decided to have some fun with rvertnet after seeing what Vijay did with Scrub Jays and Blue Jay distributions. I wanted to take a look at the distribution of HerpNET records for the classic ring-species complex Ensatina eschscholtzii. My aim was to reproduce something like this map from the excellent California Herps site:
http://www.californiaherps.com/salamanders/maps/ensatinamap.jpg
I essentially copied Vijay’s approach (see his post for details on that) but I modified to basemap to focus just on California using the “state” database in the maps package, and by manually specifying latitude and longitude limits for the map using xlim() and ylim().
I picked colors to match those used on the California Herps map. A legend would be nice, but that is left as an excercise for the reader (i.e. if you know of an easy way to do it let me know, because I was too lazy to manually set all of the legend parameters myself).
The colors used on the map correspond to: Light Blue = E. e. croceater – Yellow-blotched Ensatina; Purple = E. e. eschscholtzii – Monterey Ensatina; Blue = E. e. klauberi – Large-blotched Ensatina; Red = E. e. oregonensis – Oregon Ensatina; Black = E. e. picta – Painted Ensatina; Orange = E. e. platensis – Sierra Nevada Ensatina; E. e. xanthoptica – Yellow-eyed Ensatina.
Overall I was pretty happy with the results, though it is interesting to note that E. e. oregonensis seems to extend somewhat further south than shown on the California Herps map, overlapping completely with the western population of E. e. xanthoptica. According to California Herps, all E. e. oregonensis occuring in California are now regarded as intergrades, so this map may reflect some outdated taxonomic practices. There are also a few interesting “out-of-bounds” records, notably several E. e. eschscholtzii in Northern California far beyond the “normal” northern limit of their range near Monterey Bay. I have no idea if these might be misidentified or whether they are true vagrants.
I decided that it would be interesting to plot the location of hybrids on the map too, given that the presence/absence of hybridization zones is an important compenent of the ring species idea. It turns out that the ScienfiticName field from HerpNET was not as clean as I was hoping. When I downloaded the entire Ensatina dataset from HerpNET and checked the constitent taxa with levels() I got this:
[1] "Ensatina eschscholtzi" "Ensatina ensatina xonthoptica" [3] "Ensatina escholtzi" "Ensatina eschschlotzii eschschlotzii" [5] "Ensatina eschscholtzi" "ENSATINA ESCHSCHOLTZI" [7] "Ensatina eschscholtzi croceator" "Ensatina eschscholtzi eschscholtzi" [9] "Ensatina eschscholtzi eschscholtzi x xanthopicta" "Ensatina eschscholtzi klauberi" [11] "Ensatina eschscholtzi oregonensis" "Ensatina eschscholtzi oregonensis x xanthopicta" [13] "ENSATINA ESCHSCHOLTZI OREGONESIS" "Ensatina eschscholtzi picta" [15] "Ensatina eschscholtzi picta x oregonensis" "Ensatina eschscholtzi platensis" [17] "Ensatina eschscholtzi platensis x croceator" "Ensatina eschscholtzi xanthoptica" [19] "Ensatina eschscholtzii" "ENSATINA ESCHSCHOLTZII" [21] "Ensatina eschscholtzii cf oregonensis" "Ensatina eschscholtzii croceater" [23] "Ensatina eschscholtzii croceator" "Ensatina eschscholtzii escholtzi" [25] "Ensatina eschscholtzii eschscholtzi" "Ensatina eschscholtzii eschscholtzii" [27] "ENSATINA ESCHSCHOLTZII ESCHSCHOLTZII" "Ensatina eschscholtzii eschscholtzii x Ensatina eschscholtzii klauberi" [29] "Ensatina eschscholtzii eschscholtzii x eschscholtzii oregonensis" "Ensatina eschscholtzii eschscholtzii x xanthoptica" [31] "Ensatina eschscholtzii klauberi" "Ensatina eschscholtzii oregonensis" [33] "ENSATINA ESCHSCHOLTZII OREGONENSIS" "Ensatina eschscholtzii oregonensis x Ensatina eschscholtzii xanthoptica" [35] "Ensatina eschscholtzii oregonensis x eschscholtzii picta" "Ensatina eschscholtzii oregonensis X picta" [37] "Ensatina eschscholtzii oregonensis x platensis" "Ensatina eschscholtzii oregonensis x xanthoptica" [39] "Ensatina eschscholtzii picta" "ENSATINA ESCHSCHOLTZII PICTA" [41] "Ensatina eschscholtzii picta x oregonensis" "Ensatina eschscholtzii platensis" [43] "ENSATINA ESCHSCHOLTZII PLATENSIS" "Ensatina eschscholtzii platensis x Ensatina eschscholtzii xanthoptica" [45] "Ensatina eschscholtzii ssp." "Ensatina eschscholtzii xanthoptica" [47] "ENSATINA ESCHSCHOLTZII XANTHOPTICA" "Ensatina eschscholzii" [49] "Ensatina sp."
Yikes. Note that “Ensatina eschscholtzii oregonensis x Ensatina eschscholtzii xanthoptica” and “Ensatina eschscholtzii oregonensis x xanthoptica” appear as separate levels. This variable formatting is probably worth being aware of if you work on species that hybridize. Likewise, I noticed that my original queries also returned hybrids when the complete parent taxon name was included in taxonomic identification of the specimen. Thus an individual identified as “Ensatina eschscholtzii oregonensis x Ensatina eschscholtzii xanthoptica” would appear in queries for “Ensatina eschscholtzii oregonensis” and “Ensatina eschscholtzii xanthoptica” but a hybrid identified as “Ensatina eschscholtzii oregonensis x xanthoptica” would only appear in the former.
I began cleaning this up using gsub() to standardize the formatting when eventually my wife came in the room and asked “have you been working on that salamander thing ALL day?” At which point I realized I should probably get back to my own dissertation work.
Fun stuff though. Thanks Vijay!
library(rvertnet) library(ggplot2) library(maps) YBE<-vertoccurrence(t="Ensatina eschscholtzii croceater",grp="herp") YBE2<-subset(YBE,Latitude !=0 & Longitude != 0) ME<-vertoccurrence(t="Ensatina eschscholtzii eschscholtzii",grp="herp") ME2<-subset(ME,Latitude !=0 & Longitude != 0) LBE<-vertoccurrence(t="Ensatina eschscholtzii klauberi",grp="herp") LBE2<-subset(LBE,Latitude !=0 & Longitude != 0) OE<-vertoccurrence(t="Ensatina eschscholtzii oregonensis",grp="herp") OE2<-subset(OE,Latitude !=0 & Longitude != 0) PE<-vertoccurrence(t="Ensatina eschscholtzii picta",grp="herp") PE2<-subset(PE,Latitude !=0 & Longitude != 0) SNE<-vertoccurrence(t="Ensatina eschscholtzii platensis",grp="herp") SNE2<-subset(SNE,Latitude !=0 & Longitude != 0) YE<-vertoccurrence(t="Ensatina eschscholtzii xanthoptica",grp="herp") YE2<-subset(YE,Latitude !=0 & Longitude != 0) all_states<-map_data("state") states <- subset(all_states, region %in% c("california") ) emap <- ggplot() emap <- emap + geom_polygon( data=states, aes(x=long, y=lat, group = group),colour="white", fill="grey90" )+theme_bw() emap + geom_jitter(data = YBE2,aes(Longitude, Latitude), alpha=0.3, color = "light blue") + opts(title = "Ensatina subspecies")+ geom_jitter(data = ME2,aes(Longitude, Latitude), alpha=0.3, color = "purple")+ geom_jitter(data = LBE2, aes(Longitude, Latitude), alpha=0.3, color = "blue")+ geom_jitter(data = OE2,aes(Longitude, Latitude), alpha=0.3, color = "red")+ geom_jitter(data = PE2,aes(Longitude, Latitude), alpha=0.3, color = "black")+ geom_jitter(data = SNE2, aes(Longitude, Latitude), alpha=0.3, color = "orange")+ geom_jitter(data = YE2, aes(Longitude, Latitude), alpha=0.3, color = "yellow")+ xlim(c(-125,-113))+ylim(c(30,43))
R-bloggers.com offers daily e-mail updates about R news and tutorials about learning R and many other topics. Click here if you're looking to post or find an R/data-science job.
Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.