## Data Visualization: Basic GIS-Tasks with R  At first I want to show how you can make simple maps with R. Furthermore I want to demonstrate one important function from the package “gdalUtils”, which contains wrappers for the Geospatial Data Abstraction Library (GDAL) Utilities. I  will rasterize a vector with the function gdal_rasterize, because  it is quite fast for big vector data sets.

```library(ggplot2)
library(maptools)
library(rgdal)
library(raster)
library(plyr)
library(gdalUtils)
library(RColorBrewer)
library(rasterVis)```

Load world data (package maptools) and plot it to get an over view.

```data(wrld_simpl)
plot(wrld_simpl)``` Subset data. Select country of Austria (aut_poly) & Austria with neighbors in the west (aut_neigh).

```aut<-wrld_simpl[wrld_simpl@data\$NAME=="Austria",]
aut_neigh<-wrld_simpl[wrld_simpl@data\$SUBREGION==155,]```

Ggplot only supports data.frames, that is the reason why we transform the spatial polygon class to a data frame.

```aut_gg<-fortify(aut)
aut_neigh_gg<-fortify(aut_neigh)
p1<-ggplot()+geom_map(data=aut_neigh_gg,map=aut_neigh_gg,aes(x=long,y=lat,map_id=id),fill="white",color="black", size=0.25)+geom_map(data=aut_gg,map=aut_gg,aes(map_id=id,fill=id))
p1``` Now we change the theme of the map and add labels. Labels need a vector to define the
position. Therefore we calculate the mean X- and Y-coordinate and save both in the data.frame “labs.out_neigh”.
p2 shows us a new map with a title and labels for each country. If you want a better position for the lables you can expand the function geom_text with further arguments: Some  arguments are angle, vjust, size….

```labs.aut_neigh<-ddply(aut_neigh_gg,.(id),summarise, meanx = mean(range(long)),meany=mean(range(lat)))
p2<-p1+ggtitle("Part of Europe")+geom_text(data=labs.aut_neigh,aes(meanx,meany,label=id))
p3<-p2+theme_bw()+theme(legend.position = "right")+
scale_fill_discrete(name="Country-Id",l=50)
p3``` I am often faced with the task to transform vector data (shape files) to a raster file. There are few ways to get a solution. I have prepared two examples.

Solution one :

```extent(aut_neigh) # Delivers us the extent of our area of interest.
raster_aut_neigh_1<-raster(extent(aut_neigh)) # Creates an raster based on the extent
raster_aut_neigh_1 # Shows use the number of rows and columns. Used in the following row
raster_aut_neigh_1[] = runif(10*10) # Fills the rasters with randomly selected values.
plot(raster_aut_neigh_1) # Plot raster ```res(raster_aut_neigh_1)=c(0.5,0.5) # We increase the spatial resolution
raster_aut_neigh_1[] = runif(28*44)
plot(raster_aut_neigh_1)
``` ```raster_aut_neigh_final<-rasterize(aut_neigh,raster_aut_neigh_1,field="ISO3")
proj4string(raster_aut_neigh_final)=aut_neigh@proj4string
plot(raster_aut_neigh_final)``` Solution two:

Let’s assume that we already have a shape file. So we use the gdalUtils package and transform the shape file to a raster file. At first we create a shape file of our area of interest.

```diroutput<-"C:/Users/YourName/Documents/Foo/World_map" #define your path
# Writes a shape file into your directory
writeOGR(aut_neigh,dsn=diroutput,layer="aut_neigh",driver="ESRI Shapefile")
# Load shape file in R
```

Now we transform the shape file  (shape_file) in to a raster file.

```#find EPSG code for the Raster
EPSG <- make_EPSG() # Functions to create a data.frame which contains all EPSG-Codes
aut_neigh@proj4string@projargs # obtains the projection
# Looking after the EPSG-Codes: you will get one row with EPSG: 4326
EPSG[grepl("+proj=longlat", EPSG\$prj4, fixed=TRUE) & grepl("+datum=WGS84", EPSG\$prj4, fixed=TRUE), ]
#> 249 4326 # WGS 84 +proj=longlat +datum=WGS84 +no_defs
```

Finally we define two paths for the last calculation.  Please adapte my code to your system properties.

```inputdir<-"C:/Users/Foo/Documents/Ba/World_map"
outputdir<-"C:/Users/Foo/Documents/Ba/World_map/gdal_shp_2_r.tif"

#tr= target layer resolution; a_srs=EPSG-Code; l= shape file layer name;
#a = shape file attribute field to copy
r_id<-gdal_rasterize(inputdir, outputdir,tr=c(0.5,0.5),a_nodata=-9999,a_srs="epsg:4326",l="aut_neigh",a="UN",output_Raster = FALSE)
r_id2<-raster(outputdir)
# define raster field attribute name
names(r_id2)<-"ID"
```

Now we plot our final results.

```# We take a look to the relationship between country and the raster cell value "UN"
sort_aut<-aut_neigh@data[order(aut_neigh@data\$UN),]```

With following commands we reclassify the numeric raster call values.
Instead of the numeric cell values the legend shows us categories, which relate to the numeric attribute “UN”

```r_test <- ratify(r_id2)
rat <- levels(r_test)[]
rat\$country <- c(as.character(sort_aut[sort_aut\$UN%in%rat[,1],"ISO3"]))
levels(r_test)<-rat
my.color<-brewer.pal(7,"Dark2") # define your own color palette
levelplot(r_test,col.regions=my.color)

```  