Uncategorized

Plotting Soil Composition

#Load packages
library(GSIF)
library(soiltexture)
library(aqp)
library(plyr)
TT.plot( class.sys = "none" )
TT.plot( class.sys = "USDA.TT" )

#convert texture to classifications
vert <- TT.vertices.tbl(class.sys = "USDA.TT")
vert$x <- 1-vert$SAND+(vert$SAND-(1-vert$SILT))*0.5
vert$y <- vert$CLAY*sin(pi/3)

USDA.TT <- data.frame(TT.classes.tbl(class.sys = "USDA.TT", collapse = ", "))
TT.pnt <- as.list(rep(NA, length(USDA.TT$name)))
poly.lst <- as.list(rep(NA, length(USDA.TT$name)))
for(i in 1:length(USDA.TT$name)){
TT.pnt[[i]] <- as.integer(strsplit(unclass(paste(USDA.TT[i, "points"])), ", ")[[1]])
poly.lst[[i]] <- vert[TT.pnt[[i]],c("x","y")]
poly.lst[[i]] <- Polygons(list(Polygon(rbind(poly.lst[[i]], poly.lst[[i]][1,]))), ID=i)
}

#convert texture triangle to a spatial object:
poly.sp <- SpatialPolygons(poly.lst, proj4string=CRS(as.character(NA)))
poly.USDA.TT <- SpatialPolygonsDataFrame(poly.sp, data.frame(ID=USDA.TT$name), match.ID=FALSE)
spplot(poly.USDA.TT)
ttcolor
#Load soil composition data
soil <- read.csv("soil.csv")
TT.plot(
class.sys = "HYPRES.TT",
tri.data = soil,
main = "Soil texture data"
)
tt
TT.plot(class.sys = "none",
tri.data = soil,
z.name = "Organic Content",
main = "Soil texture triangle and Organic Content bubble plot"
)

tt2








Maps

R foR Rasters

require(rgdal)
require(raster)
require(ggplot2)
require(dplyr)
require(rasterVis)

setwd("C:/DATA/Blog/raster/NLCD2011_LC_N39W072")
li <- raster("NLCD2011_LC_N39W072.tif")
plot(li)
li
shape <- readOGR(dsn = "C:/DATA/Blog/raster/lishape/government_units_NRCSCNTY_ny_3451580_01/government_units", layer = "county_nrcs_a_ny")
suffolk <- subset(shape, COUNTYNAME == "Suffolk")
plot(suffolk)
suf
suffolk_c <- crop(li, extent(suffolk))
suffolk_c2 <- mask(suffolk_c, suffolk)
suff

dist <- data.frame(freq(suffolk_c2))
g <- ggplot(dist, aes(x=X.1, y=count2, fill=X.1)) + geom_bar(stat = "identity") +
 coord_flip()

hist2
Uncategorized

Word Clouds

Word clouds are a cool way to show themes in text.  I did this word cloud with my master's thesis!

1.Require packages
 library(RXKCD)
 library(tm)
 library(wordcloud)
 library(RColorBrewer)
 library(SnowballC)

2.Pull out text file
 filePath <- "C:/data/blog/wordcloud/kk.txt"
 text <- readLines(filePath)
 docs <- Corpus(VectorSource(text))
 inspect(docs)

3.Get rid of messy characters
 toSpace <- content_transformer(function (x , pattern ) gsub(pattern, " ", x))
 docs <- tm_map(docs, toSpace, "/")
 docs <- tm_map(docs, toSpace, "@")
 docs <- tm_map(docs, toSpace, "\\|")

4.Convert the text to lower case
 docs <- tm_map(docs, content_transformer(tolower))

5. Remove numbers
 docs <- tm_map(docs, removeNumbers)

6. Remove english common stopwords
 docs <- tm_map(docs, removeWords, stopwords("english"))

7. Remove punctuations
 docs <- tm_map(docs, removePunctuation)

8.Eliminate extra white spaces
 docs <- tm_map(docs, stripWhitespace)

9. Convert to matrix
 dtm <- TermDocumentMatrix(docs)
 m <- as.matrix(dtm)
 v <- sort(rowSums(m),decreasing=TRUE)
 d <- data.frame(word = names(v),freq=v)
 head(d, 10)

10. Create wordcloud
 wordcloud(words = d$word, freq = d$freq, min.freq = 1,
 max.words=200, random.order=FALSE, rot.per=0.35,
 colors=brewer.pal(8, "Dark2"))

pic

Maps

Coral Reef Bleaching

1.Start with some coral reef data
 library(dplyr)
 library(ggplot2)
 library (rgdal)
 library (rgeos)
 library(maptools)
 library(tmap)

cb <- read.csv("CoralBleaching.csv")

2.Read in reef data with latitudes and longitudes
 b <- ggplot(cb, aes(LON, LAT)) + coord_cartesian(ylim=c(-40, 40), xlim=c(-200, 200))
 b + geom_point()
 b + geom_point(aes(color=cb$BLEACHING_SEVERITY))+scale_colour_brewer("Bleaching Severity", palette="PiYG")+ coord_equal(ratio=1)

reefvleach1

3. Add in a shapefile of coastlines
 coast <-readOGR("ne_50m_coastline.shp", layer="ne_50m_coastline")
 ggplot() + geom_polygon(data=coast, aes(x=long, y=lat, group=group))

4. Add your point data on the map and color according to bleaching severity
 ggplot() + geom_polygon(data=coast aes(x=long, y=lat, group=group), fill="grey40",
 colour="grey90", alpha=1)+labs(x="", y="", title="Coral Reef Bleaching")+ #labels
 theme(axis.ticks.y = element_blank(),axis.text.y = element_blank(),
 plot.title = element_text(lineheight=.8, face="bold", vjust=1))+space
 geom_point(aes(x=LON, y=LAT, color=BLEACHING_SEVERITY), data=cb, alpha=1, size=3, color="grey20")+geom_point(aes(x=LON, y=LAT, color=BLEACHING_SEVERITY), data=cb, alpha=1, size=2)+scale_colour_brewer("Bleaching Severity",
 palette="PiYG")+ coord_equal(ratio=1)

reefbleachpic

http://zevross.com/blog/2014/07/16/mapping-in-r-using-the-ggplot2-package/
http://www.sthda.com/english/wiki/ggplot2-colors-how-to-change-colors-automatically-and-manually
http://remi-daigle.github.io/GIS_mapping_in_R/

Maps

Google Maps & R Maps

You can grab satellite imagery from google maps by simply calling up the area by latitude and longitude

  1. Load packages
    library(ggmap)
    library(ggplot2)
  2. Set a range
    lat <- c(18, 23)
    lon <- c(-161, -154)
  3. Get a map
    map <- get_map(location = c(lon = mean(lon), lat = mean(lat)), zoom = 6, maptype = “satellite”, source = “google”)
    plot(map)
    hi
    4.  You can also play with the zoom and image type ie :zoom into Oahulat <- c(21.1, 21.8)
    lon <- c(-158.3, -157.567)

    map <- get_map(location = c(lon = mean(lon), lat = mean(lat)), zoom = 10, maptype = “satellite”, source = “google”)
    plot(map)
    oahu
    map2 <- get_map(location = c(lon = mean(lon), lat = mean(lat)), zoom = 10, maptype = “terrain”, source = “google”)
    plot(map2)
    oahu_terrain
    and if you want to make your map really pretty…
    map3 <- get_map(location = c(lon = mean(lon), lat = mean(lat)), zoom = 10, maptype = “watercolor”, source = “google”)
    plot(map3)

    oahu_watercolor
    https://www.r-bloggers.com/google-maps-and-ggmap/

Uncategorized

Census Bureau Maps

Maps

1. Download population data from Census Bureau
I chose population data by county for 2015 American Community Survey
download and save (aff_download.zip)
2. Download shapefile of counties in United States (tl_2016_us_county.zip)
3. Extract compressed files

4. Load required packages
require(dplyr)
require(geosphere)
require(ggmap)
require(ggplot2)
require(maptools)
require(raster)
require(rgdal)
require(rgeos)
require(sp)
require(spdep)

5. Establish directories
xdir = ‘C:/_data/blog/map/Shapefiles’
ydir = ‘C:/_data/blog/map’

6. Read spatial data into R by directory and filename (object is currently a spatial polygon dataframe)
SPD <- readOGR(dsn = xdir, layer = ‘tl_2016_us_county’)
SPD@data$id <- rownames(SPD@data)

7. Covert to dataframe
DF <- fortify(SPD, region= ‘id’)
DF2 <- merge(DF, SPD@data, by=’id’)

8. Filter by desired latitude and longitude
DF3 <- filter(DF2, lat < 48)
DF4 <- filter(DF3, lat > 34)
DF5 <- filter(DF4, long > -76)
DF6 <- filter(DF5, long < -67)

9. Change working directory to where population data is stored and read first 4 columns
setwd(ydir)
popdata <- read.csv(“ACS_15_5YR_S0101_with_ann.csv”)
popdata2 <- popdata[,1:4]

10. Join spatial data and population data
map <- left_join(DF6, popdata2, by = c(“GEOID” = “GEO.id2”))

11. Convert population to numeric so scale will be continuous
map$pop <- as.numeric(map$HC01_EST_VC01)

12. Plot map dataframe using ggplot with latitude and logitude being x and y
p <- ggplot() + geom_polygon(data=test, aes(x=long, y=lat, group=group, fill=pop))
p2 <- p + coord_equal()
print(p2)

map1
13. You can also change the map colors to better visualize the data
Here I made counties with high populations red and low populations green.  We could use this type of visualization if we are discussing a simple concept such as high populations having negative impacts on ecosystem health

map2
References:
http://ggplot2.tidyverse.org/reference/scale_brewer.html
http://www.cookbook-r.com/Graphs/Colors_(ggplot2)/

summarystatistics

Summary Statistics

The benefit of working with R as opposed to excel is that you can work easily with large amounts of data.  The first thing you want to do with a large data set is to observe the structure of the data as well as summary statistics.

1.Grab data
I grabbed some data from world bank about threatened bird species

1.Set working directory
setwd(“C:/Users/kkuntz/_data/blog/bird”)

2.Read data
dat <- read.csv(“birds.csv”)

3.Filter for threatened bird species
dat2 <- filter(dat, Indicator_Name==”Bird species, threatened”)

4. Take a look at your data structure
nrow(dat2)
str(dat2)

5. Perform some summary statistics
var(dat2$Z2016) = variance
sd(dat2$Z2016) = standard deviation
max(dat2$Z2016) = maximum
min(dat2$Z2016) = minimum
range(dat2$Z2016) = range

6. Plot a histogram to visualize data
I filtered by 200 to the areas with the highest threatened species

dat5 <- filter(dat4, X2016>200)

g <- ggplot(data=dat5, aes(x=Country.Name, y=X2016, fill = Country.Name))+ geom_bar(stat=”identity”) +guides(fill=FALSE)

g + theme(axis.text.y = element_text(size=7, angle = 0, hjust=0)) + coord_flip()

g + theme(axis.text.x = element_text(size=15, angle = 90, hjust=0)+ theme(axis.text.y = element_text(size=15, angle = 90, hjust=0)) + coord_flip())

graph3