data(Test_Data_US)
g <- GeoDistPSU(Test_Data_US$lat,
Test_Data_US$long,
"miles", 100,
Input.ID = Test_Data_US$ID)
# Plot GeoDistPSU output
plot(g$PSU.Info$PSU.Mean.Longitude,
g$PSU.Info$PSU.Mean.Latitude,
col = 1:nrow(g$PSU.Info),
pch = 19,
main = "Plot of PSU Centers",
xlab = "Longitude",
ylab = "Latitude")
grid(col = "grey40")
# Plot GeoDistPSU output with map
if (FALSE) {
# install package sf to run usmap_transform
library(ggplot2)
library(sp)
library(usmap)
# Transform PSUs into usmap projection
g.map <- cbind(long = g$PSU.Info$PSU.Mean.Longitude,
lat = g$PSU.Info$PSU.Mean.Latitude)
g.map <- as.data.frame(g.map)
g.proj <- usmap::usmap_transform(g.map,
input_names = c("long", "lat"),
output_names = c("Long", "Lat"))
usmap::plot_usmap(color = "gray") +
geom_point(data = g.proj,
aes(x = Long,
y = Lat))
# Create histogram of maximum distance
hist(g$PSU.Info$PSU.Max.Dist,
main = "Histogram of Maximum Within-PSU Distance",
xlab = "Distance",
ylab = "Frequency")
}
Run the code above in your browser using DataLab