R - Applying the same script to multiple dataframes

Question

I have a script which reads a .dat file containing co-ordinates, rainfall and crop area data as below.

  East  North    rain  Wheat Sbarley Potato OSR fMaize Total  LCA
10000 510000 1498.73 0.0021     0.5 0.0022   0 0.0056  0.01 0.01
10000 510000 1498.73 0.0021   0.034 0.0022   0 0.0056  0.01 0.01
10000 510000 1498.73 0.0021   0.001 0.0022   0 0.0056  0.01 0.01
10000 515000 1518.51 0.0000    0.12 0.0000   0 0.0000  0.00 0.00
10000 515000 1518.51 0.0000  0.0078 0.0125   0 0.0000  0.00 0.00
10000 515000 1518.51 0.0000       0 0.0000   0   0.03  0.00 0.00

The code below calculates the emissions for Wheat through a series of models before extracting the data and producing a raster file and plotting it in ggplot. Having read these related queries I'm getting rather confused as to whether I need to make a package or am missing something very basic as to how to repeat the code through each crop type.

library(doBy)
library(ggplot2)
library(plyr)
library(raster)
library(RColorBrewer)
library(rgdal)
library(scales)
library(sp)

### Read in the data
crmet <- read.csv("data.dat")
# Remove NA values
crm <- crmet[ ! crmet$rain %in% -119988, ]
crm <- crm[ ! crm$Wheat %in% -9999, ]

### Set model parameters
a <- 0.1474
b <- 0.0005232
g <- -0.00001518
d <- 0.000003662
N <- 182

### Models
crm$logN2O <- a+(b*crm$rain)+(g*N)+(d*crm$rain*N)
crm$eN2O <- exp(crm$logN2O)
crm$whN2O <- crm$eN2O*crm$Wheat

### Prepare data for conversion to raster
crmet.ras <- crm
crmet.ras <- rename(crmet.ras, c("East"="x","North"="y"))

#### Make wheat emissions raster
wn <- crmet.ras[,c(1,2,13)]
spg <- wn

# Set the Eastings and Northings to coordinates
coordinates(spg) <- ~ x + y
# coerce to SpatialPixelsDataFrame
gridded(spg) <- TRUE
# coerce to raster
rasterDF <- raster(spg)
# Add projection to it - in this case OSBG36
proj4string(rasterDF) <- CRS("+proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 +y_0=-100000 +ellps=airy +datum=OSGB36 +units=m +no_defs")
rasterDF

writeRaster(rasterDF, 'wn.tif', overwrite=T)

### Plotting the raster:
whplot <-ggplot(wn, aes(x = x, y = y))+
  geom_tile(aes(fill = whN2O))+
  theme_minimal()+
  theme(plot.title = element_text(size=20, face="bold"),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        axis.ticks = element_blank(),
        axis.text.x = element_blank(),
        axis.text.y = element_blank(),
        legend.title = element_text(size=16, face="bold"))+
  scale_fill_gradient(name=expression(paste("", N[2], "O ", Ha^-1, sep="")))+
  xlab ("")+
  ylab ("")+
  labs (title="Nitrous oxide emissions\nfrom Wheat")

whplot

I'm trying to turn as much of the above into functions as I can but they're all rather more complicated than the examples I find on here and in the ?help files. Many thanks in advance for any help/suggestions.

citronrose · Accepted Answer · 2016-08-18 12:49:19Z

1

You need to put your code in a function that takes the data file name as argument. Then you can call your function with the relevant data file name. Something like:

library(doBy)
library(ggplot2)
library(plyr)
library(raster)
library(RColorBrewer)
library(rgdal)
library(scales)
library(sp)

#defining the function
my.neat.function <- function(datafname){
### Read in the data
crmet <- read.csv(datafname)
# Remove NA values
crm <- crmet[ ! crmet$rain %in% -119988, ]
crm <- crm[ ! crm$Wheat %in% -9999, ]

### Set model parameters
a <- 0.1474
b <- 0.0005232
g <- -0.00001518
d <- 0.000003662
N <- 182
### Models
crm$logN2O <- a+(b*crm$rain)+(g*N)+(d*crm$rain*N)
crm$eN2O <- exp(crm$logN2O)
crm$whN2O <- crm$eN2O*crm$Wheat

### Prepare data for conversion to raster
crmet.ras <- crm
crmet.ras <- rename(crmet.ras, c("East"="x","North"="y"))

#### Make wheat emissions raster
wn <- crmet.ras[,c(1,2,13)]
spg <- wn

# Set the Eastings and Northings to coordinates
coordinates(spg) <- ~ x + y
# coerce to SpatialPixelsDataFrame
gridded(spg) <- TRUE
# coerce to raster
rasterDF <- raster(spg)
# Add projection to it - in this case OSBG36
proj4string(rasterDF) <- CRS("+proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 +y_0=-100000 +ellps=airy +datum=OSGB36 +units=m +no_defs")
rasterDF

writeRaster(rasterDF, 'wn.tif', overwrite=T)

### Plotting the raster:
whplot <-ggplot(wn, aes(x = x, y = y))+
  geom_tile(aes(fill = whN2O))+
  theme_minimal()+
  theme(plot.title = element_text(size=20, face="bold"),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        axis.ticks = element_blank(),
        axis.text.x = element_blank(),
        axis.text.y = element_blank(),
        legend.title = element_text(size=16, face="bold"))+
  scale_fill_gradient(name=expression(paste("", N[2], "O ", Ha^-1, sep="")))+
  xlab ("")+
  ylab ("")+
  labs (title="Nitrous oxide emissions\nfrom Wheat")

whplot
} #end of function definition
my.neat.function("data.dat") #first call to function
my.neat.function("otherdata.dat")#same thing with another dataset

If the model parameters for different data, you need to add the vector of parameter values as argument to the function.

answered Aug 18, 2016 at 12:49

citronrose

4093 silver badges9 bronze badges

Sign up to request clarification or add additional context in comments.

2 Comments

B.Wel Over a year ago

Cool, thanks! If I subset the data into one for each species (Wheat, SBarley, etc.) if I use crm[[3]] in place of crm$Wheat will that read the third column in each subset dataframe?

B.Wel Over a year ago

Realized the easiest fix to get the function to work was to simplify the original data frame to just five columns using the melt function from reshape2.

Collectives™ on Stack Overflow

R - Applying the same script to multiple dataframes

1 Answer 1

2 Comments

Your Answer

Linked

Hot Network Questions

Collectives™ on Stack Overflow

1 Answer 1

2 Comments

Your Answer

Sign up or log in

Post as a guest

Linked

Related