Skip to content

Instantly share code, notes, and snippets.

@johnlees
Forked from epijim/jb_pubdistance.r
Last active August 29, 2015 14:11
Show Gist options
  • Save johnlees/e62f95d7f96432412a3f to your computer and use it in GitHub Desktop.
Save johnlees/e62f95d7f96432412a3f to your computer and use it in GitHub Desktop.
jb_pubdistance <- function(
v_pubs=c("The Cambridge Brew House","The Cambridge Blue",
"King Street Run P.H.","The Elm Tree","The Regal"),
v_location="Cambridge, UK",
v_zoom=13,
listpubs=FALSE,
cam_pubs=TRUE,
crow_distances=FALSE,
units="minutes"){
## Read in pub data
temporaryFile <- tempfile()
download.file("https://gist.githubusercontent.com/epijim/47da5037b6bb0d750b9b/raw/1f5608b826d039384775f2de87fcab9f5c8311d9/cam_pubs.csv",
destfile=temporaryFile, method="curl")
pubs <- read.csv(temporaryFile)
if(listpubs==TRUE){
print(pubs$name)
stop("listpubs TRUE: function stopped and pubs listed")
}
# required package
library(ggmap)
library(hash)
# load map
basemap <- get_map(location=v_location,
zoom=v_zoom, source='google',
maptype="terrain", color="bw")
if(cam_pubs==TRUE){
if(crow_distances==FALSE){
## Read in distance data
temporaryFile <- tempfile()
download.file("https://gist.githubusercontent.com/epijim/1b049708cef27d9ddc34/raw/02a4f8559541b7dd5589ef1f530966759d81a711/campubdistances.csv",
destfile=temporaryFile, method="curl")
googled_distances <- read.csv(temporaryFile)
distance_hash <- hash(keys=paste(googled_distances$latlon_1,googled_distances$latlon_2),
values=googled_distances$minutes)
pubs$latlon <- paste0(pubs$lat,",",pubs$lon)
## Choose pubs you want #########################################
selected.pubs <- subset(pubs, name %in% v_pubs)
## Variables to set #############################################
n <- as.numeric(nrow(selected.pubs)) # number of pubs
M <- 5000 # number of iterations
temp <- 100 # initial temperature
finaltemp <- 0.1
tempfactor <- (temp/finaltemp)^(-1/M)
fixed.itr <- 20 # number of iteration for a fixed temperature
## FUNCTIONS ####################################################
## Check if integer
check.integer <- function(x) {
x == round(x)
}
## Difference between 2 points
jb_difference <- function(testv1,testv2){
d <- distance_hash[[paste(testv1,testv2)]]
if(is.null(d))
{
d <- distance_hash[[paste(testv2,testv1)]]
}
d <- as.numeric(d)
return(d)
}
distance <- function(ord){
p <- length(ord)
d <- rep(0,p)
for(i in 1:(p-1)){
d[i] <- jb_difference(selected.pubs$latlon[ord[i]],
selected.pubs$latlon[ord[i+1]])
}
d[p] <- jb_difference(selected.pubs$latlon[ord[1]],
selected.pubs$latlon[ord[p]] )
distance=sum(d)
return(distance)
}
## Initial values ################################################
cost1 <- rep(0,M)
numaccept <- 0 # to keep track of accepting
templist <- rep(0,M) # to keep track of temperature
ord <- seq(1:n) # initial order for the path
max.distance <- 15
## Run the loop ################################################
for(i in 1:M){
cost1[i]=distance(ord)
if (i %% 10 == 0)
{
cat(sprintf("Iteration %d of %d\n",i,M))
}
for(t in 1:(fixed.itr)){
# propose to move
a=sample(1:n,2) # randomly choosing two cities
result=ord
result[a[2]]=ord[a[1]] # swapping
result[a[1]]=ord[a[2]] # swapping
U=runif(1)
if(
U < exp( (distance(ord)-distance(result))/temp )){
ord<-result # accept/reject
numaccept <- numaccept + 1 # to compute the acceptance rate +}
}
}
templist[i] <- temp
temp <- temp*0.999
temp <- tempfactor*temp
}
temperature <- data.frame(1:M,templist)
length <- 1:length(cost1)
fitting <- data.frame(length,cost1)
pubs_inorder <- selected.pubs[ord,]
map <- ggmap(basemap) +
geom_segment(aes(xend=c(tail(lon, n=-1), NA),
yend=c(tail(lat, n=-1), NA)),
data=pubs_inorder)
result <- list(distance=cost1[M],
pubs_inorder=pubs_inorder,
fittingprocess=fitting,
temperature=temperature)
}
if(crow_distances==TRUE){
## Choose pubs you want #########################################
selected.pubs <- subset(pubs, name %in% v_pubs)
## Variables to set #############################################
n <- as.numeric(nrow(selected.pubs)) # number of pubs
M <- 4000 # number of iterations
temp <- 100 # initial temperature
finaltemp <- 0.1
tempfactor <- (temp/finaltemp)^(-1/M)
fixed.itr <- 20 # number of iteration for a fixed temperature
## FUNCTIONS ####################################################
## Check if integer
check.integer <- function(x) {
x == round(x)
}
## Difference between 2 points
jb_difference <- function(lat_1,lon_1,lat_2,lon_2){
rad <- pi/180
a1 <- lat_1*rad
a2 <- lon_1*rad
b1 <- lat_2*rad
b2 <- lon_2*rad
dlon <- b2 - a2
dlat <- b1 - a1
a <- (sin(dlat/2))^2 + cos(a1) * cos(b1) * (sin(dlon/2))^2
c <- 2 * atan2(sqrt(a), sqrt(1 - a))
R <- 6378.14500001
d <- R * c
return(d)
}
distance <- function(ord){
p <- length(ord)
d <- rep(0,p)
for(i in 1:(p-1)){
d[i] <- jb_difference(selected.pubs$lat[ord[i]],
selected.pubs$lon[ord[i]],
selected.pubs$lat[ord[i+1]],
selected.pubs$lon[ord[i+1]])
}
d[p] <- jb_difference(selected.pubs$lat[ord[1]],
selected.pubs$lon[ord[1]],
selected.pubs$lat[ord[p]],
selected.pubs$lon[ord[p]] )
distance=sum(d)
return(distance)
}
## Initial values ################################################
cost1 <- rep(0,M)
numaccept <- 0 # to keep track of accepting
templist <- rep(0,M) # to keep track of temperature
ord <- seq(1:n) # initial order for the path
max.distance <- 15
## Run the loop ################################################
for(i in 1:M){
cost1[i]=distance(ord)
for(t in 1:(fixed.itr)){
# propose to move
a=sample(1:n,2) # randomly choosing two cities
result=ord
result[a[2]]=ord[a[1]] # swapping
result[a[1]]=ord[a[2]] # swapping
U=runif(1)
if(
U < exp( (distance(ord)-distance(result))/temp )){
ord<-result # accept/reject
numaccept <- numaccept + 1 # to compute the acceptance rate +}
}
}
templist[i] <- temp
temp <- temp*0.999
temp <- tempfactor*temp
}
temperature <- data.frame(1:M,templist)
length <- 1:length(cost1)
fitting <- data.frame(length,cost1)
pubs_inorder <- selected.pubs[ord,]
map <- ggmap(basemap) +
geom_segment(aes(xend=c(tail(lon, n=-1), NA),
yend=c(tail(lat, n=-1), NA)),
data=pubs_inorder)
result <- list(distance=cost1[M],
pubs_inorder=pubs_inorder,
fittingprocess=fitting,
temperature=temperature)
}
}
if(cam_pubs==FALSE){
## Variables to set #############################################
n <- as.numeric(nrow(v_pubs)) # number of pubs
M <- 4000 # number of iterations
temp <- 100 # initial temperature
finaltemp <- 0.1
tempfactor <- (temp/finaltemp)^(-1/M)
fixed.itr <- 20 # number of iteration for a fixed temperature
## FUNCTIONS ####################################################
## Check if integer
check.integer <- function(x) {
x == round(x)
}
## Difference between 2 points
jb_difference <- function(lat_1,lon_1,lat_2,lon_2){
rad <- pi/180
a1 <- lat_1*rad
a2 <- lon_1*rad
b1 <- lat_2*rad
b2 <- lon_2*rad
dlon <- b2 - a2
dlat <- b1 - a1
a <- (sin(dlat/2))^2 + cos(a1) * cos(b1) * (sin(dlon/2))^2
c <- 2 * atan2(sqrt(a), sqrt(1 - a))
R <- 6378.14500001
d <- R * c
return(d)
}
distance <- function(ord){
p <- length(ord)
d <- rep(0,p)
for(i in 1:(p-1)){
d[i] <- jb_difference(v_pubs[ord[i],1],
v_pubs[ord[i],2],
v_pubs[ord[i+1],1],
v_pubs[ord[i+1],2])
}
d[p] <- jb_difference(v_pubs[ord[1],1],
v_pubs[ord[1],2],
v_pubs[ord[p],1],
v_pubs[ord[p],2] )
distance=sum(d)
return(distance)
}
## Initial values ################################################
cost1 <- rep(0,M)
numaccept <- 0 # to keep track of accepting
templist <- rep(0,M) # to keep track of temperature
ord <- seq(1:n) # initial order for the path
max.distance <- 15
## Run the loop ################################################
for(i in 1:M){
cost1[i]=distance(ord)
for(t in 1:(fixed.itr)){
# propose to move
a=sample(1:n,2) # randomly choosing two cities
result=ord
result[a[2]]=ord[a[1]] # swapping
result[a[1]]=ord[a[2]] # swapping
U=runif(1)
if(
U < exp( (distance(ord)-distance(result))/temp )){
ord<-result # accept/reject
numaccept <- numaccept + 1 # to compute the acceptance rate +}
}
}
templist[i] <- temp
temp <- temp*0.999
temp <- tempfactor*temp
}
temperature <- data.frame(1:M,templist)
length <- 1:length(cost1)
fitting <- data.frame(length,cost1)
pubs_inorder <- v_pubs[ord,]
map <- ggmap(basemap) +
geom_segment(aes(xend=c(tail(lon, n=-1), NA),
yend=c(tail(lat, n=-1), NA)),
data=pubs_inorder)+
xlab("Latitude")+
ylab("Longitude")
result <- list(distance=cost1[M],
pubs_inorder=pubs_inorder,
fittingprocess=fitting,
temperature=temperature)
}
print(map)
return(result)
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment