I'm trying to plot some path/connection based maps but am unable to figure out how.
I see a lot of possibilities for one point based metrics (crime hotspots in London, etc. with googleVis, ggmap, etc.) but I can't find too many examples of two points based metrics (immigrations between cities, train routes, etc.) There is an example with the package geosphere
, but it seems to not be available for R 3.0.2
.
Ideally I would like something like this D3 example and would also like to customise the thickness, colour, etc. of the lines and the circles.
PS: I don't suppose rCharts
can do this just yet, right?
Here's something I started last year and never got around to polishing properly but hopefully it should answer your question of joining points on a map with great circle lines, with the flexibility to customise lines, circles etc.
I use the (my) rworldmap package for mapping, WDI for world bank data and geosphere for the great circle lines. The aim was to plot aid flows from all donor countries to all recipient countries (one plot per donor). Below is an example plot, and below that the code. Hope it helps. Would be nice to find time to pick it up again!
Andy
library(rworldmap)
library(WDI) # WORLD BANK INDICATORS
## lines of either type may obscure more than they add
##**choose line option here
addLines <- 'gc' #'none''straight' 'gc'
if ( addLines == 'gc' ) library(geosphere)
# setting background colours
oceanCol = rgb(7,0,30,maxColorValue=255)
landCol = oceanCol
#produces a list of indicator IDs and names as a matrix
indicatorList <- WDIsearch('aid flows')
#setting up a world map shaped plot window
#*beware this is windows specific
mapDevice('windows',width=10,height=4.5)
year <- 2000
#for(indNum in 1:2)
for(indNum in 1:nrow(indicatorList))
{
indID <- indicatorList[indNum][1]
donorISO3 <- substr(indID,start=8,stop=10)
dFdonor <- WDI(indicator=indID,start=year,end=year)
#divide by 10^6 for million dollars
dFdonor[indID] <- dFdonor[indID] * 1/10^6
sPDFdonor <- joinCountryData2Map(dFdonor,nameJoinColumn='country',joinCode='NAME')
#take out Antarctica
sPDFdonor <- sPDFdonor[-which(row.names(sPDFdonor)=='Antarctica'),]
legendTitle=paste("aid flow from",donorISO3,year,"(millions US$)")
mapBubbles(sPDFdonor, nameZSize=indID, plotZeroVals=FALSE, legendHoriz=TRUE, legendPos="bottom", fill=FALSE, legendTitle=legendTitle, oceanCol=oceanCol, landCol=landCol,borderCol=rgb(50,50,50,maxColorValue=255),lwd=0.5,lwdSymbols=1)
#removed because not working , main=paste('donor', donorISO3,year)
#now can I plot lines from the centroid of the donor to the centroids of the recipients
xDonor <- sPDFdonor$LON[ which(sPDFdonor$ISO3==donorISO3) ]
yDonor <- sPDFdonor$LAT[ which(sPDFdonor$ISO3==donorISO3) ]
xRecips <- sPDFdonor$LON[ which(sPDFdonor[[indID]] > 0) ]
yRecips <- sPDFdonor$LAT[ which(sPDFdonor[[indID]] > 0) ]
amountRecips <- sPDFdonor[[indID]][ which(sPDFdonor[[indID]] > 0) ]
## straight lines
if ( addLines == 'straight' )
{
for(line in 1:length(xRecips))
{
#col <- 'blue'
#i could modify the colour of the lines by the size of the donation
#col=rgb(1,1,1,alpha=amountRecips[line]/max(amountRecips))
#moving up lower values
col=rgb(1,1,0,alpha=sqrt(amountRecips[line])/sqrt(max(amountRecips)))
lines(x=c(xDonor,xRecips[line]),y=c(yDonor,yRecips[line]),col=col, lty="dotted", lwd=0.5) #lty = "dashed", "dotted", "dotdash", "longdash", lwd some devices support <1
}
}
## great circle lines
## don't work well when donor not centred in the map
## also the loop fails at CEC & TOT because not ISO3 codes
if ( addLines == 'gc' & donorISO3 != "CEC" & donorISO3 != "TOT" )
{
for(line in 1:length(xRecips))
{
#gC <- gcIntermediate(c(xDonor,yDonor),c(xRecips[line],yRecips[line]), n=50, breakAtDateLine=TRUE)
#30/10/13 lines command failed with Error in xy.coords(x, y) :
#'x' is a list, but does not have components 'x' and 'y'
#adding sp=TRUE solved
gC <- gcIntermediate(c(xDonor,yDonor),c(xRecips[line],yRecips[line]), n=50, breakAtDateLine=TRUE, sp=TRUE)
#i could modify the colour of the lines by the size of the donation
#col=rgb(1,1,1,alpha=amountRecips[line]/max(amountRecips))
#moving up lower values
col=rgb(1,1,0,alpha=sqrt(amountRecips[line])/sqrt(max(amountRecips)))
lines(gC,col=col,lwd=0.5)
}
}
#adding coasts in blue looks nice but may distract
data(coastsCoarse)
plot(coastsCoarse,add=TRUE,col='blue')
#repeating mapBubbles with add=T to go on top of the lines
mapBubbles(sPDFdonor, nameZSize=indID, plotZeroVals=FALSE, fill=FALSE, addLegend=FALSE, add=TRUE, ,lwd=2)
#removed because not working : , main=paste('donor', donorISO3,year)
#looking at adding country labels
text(xRecips,yRecips,sPDFdonor$NAME[ which(sPDFdonor[[indID]] > 0) ],col=rgb(1,1,1,alpha=0.3),cex=0.6,pos=4) #pos=4 right (1=b,2=l,3=ab)
#add a title
nameDonor <- sPDFdonor$NAME[ which(sPDFdonor$ISO3==donorISO3) ]
mtext(paste("Aid flow from",nameDonor,year), cex = 1.8, line=-0.8)
#savePlot(paste("C:\\rProjects\\aidVisCompetition2012\\Rplots\\greatCircles\\wdiAidFlowLinesDonor",donorISO3,year,sep=''),type='png')
#savePlot(paste("C:\\rProjects\\aidVisCompetition2012\\Rplots\\greatCircles\\wdiAidFlowLinesDonor",donorISO3,year,sep=''),type='pdf')
} #end of indNum loop