R Markdown for Duke Farms - Grassland Field (skeet shoot field)

Raw Data: CLICK HERE TO DOWNLOAD the raw data.

library(gstat)




 

Part 1: When ‘n = 150’


 

1.1 Compute undirected (isotropic) empirical semivariogram
var1 <- variogram(Total_C ~ 1, ~x.axis + y.axis, cutoff = 500, width = 30, data = mc)

# create figure of empirical semivariogram
plot(var1, plot.numbers = TRUE, xlab = "Distance", ylab = "Semivariance", cex = 1, cex.axis = 2, xlim=c(0, 550), ylim=c(0.0, 0.25))

Figure 1: Undirected semivariogram of Total_C .


 

1.1.2 Semivariogram models

To fit theoretical semivariograms, the function ‘fit.variogram’ of the package ‘gstat’ is used. Sill, nugget and range are set to be calculated based on the empirical variogram data, which is also used to fit the model ‘(fit.method=1)’. This method fits the variogram model to the experimental variogram, using weighted least squares with weight = Nj, where Nj is the number of observations in the j -th distance class (bin) (from: http://www.gstat. org/gstat.pdf, Table 4.2). The exponential model is used here. Semivariogram models are only fitted to undirected empirical semivariograms as the number of point pairs per bin in the directed ones is very low and predictions therefore have less power (the number of points per bin (np) for the undirected semivariograms is as high as for all directed semivariograms together).
 

compute semivariogram model based on the empiral semivariogram ‘var1’**
mod1 <- fit.variogram(var1, vgm(psill = NA, "Exp", range = NA, 1), fit.sills = TRUE, fit.ranges = TRUE, fit.method = 1)
## Warning in fit.variogram(var1, vgm(psill = NA, "Exp", range = NA, 1), fit.sills
## = TRUE, : No convergence after 200 iterations: try different initial values?

 

show nugget, psill and range of the semivariogram model
mod1
##   model      psill    range
## 1   Nug 0.04865517  0.00000
## 2   Exp 0.08499639 40.38164

 

show the (weighted) sum of squared errors of the fitted model
attr(mod1, "SSErr")
## [1] 4.346393


 

plot figure of semivariogram model
plot(var1, plot.numbers = TRUE, model = mod1, xlab = "Distance", ylab = "Semivariance", xlim=c(0, 550), ylim=c(0.0, 0.25))

Figure 2: Exponential model fit on Total_C.




 

Part 2: when ‘n = 94’ (I removed the data from 0.5-m sampling nest)

 

# select all data where the distance is not (! = is not) equal to 0.5

mc1 = subset(mc, !distance== "0.5")



2.1 Compute undirected (isotropic) empirical semivariogram
var1 <- variogram(Total_C ~ 1, ~x.axis + y.axis, cutoff = 500, width = 30, data = mc1)
#
#
# create figure of empirical semivariogram
plot(var1, plot.numbers = TRUE, xlab = "Distance", ylab = "Semivariance", cex = 1, cex.axis = 2, xlim=c(0, 550), ylim=c(0.0, 0.25))

Figure 3: Undirected semivariogram of Total_C . The empirical semivariogram shows a spatial distribution without the 0.5-m samples.



2.1.2 Semivariogram models

Now we fit theoretical semivariograms to those empirical semivariograms which showed an autocorrelation structure in the previous data analyses.


compute semivariogram model based on the empiral semivariogram ‘var1’
mod1 <- fit.variogram(var1, vgm(psill = NA, "Exp", range = NA, 1), fit.sills = TRUE, fit.ranges = TRUE, fit.method = 1)


show nugget, psill and range of the semivariogram model
mod1
##   model      psill    range
## 1   Nug 0.01721994  0.00000
## 2   Exp 0.11024465 28.24095


show the (weighted) sum of squared errors of the fitted model
attr(mod1, "SSErr")
## [1] 1.364782



plot figure of semivariogram model
plot(var1, plot.numbers = TRUE, model = mod1, xlab = "Distance", ylab = "Semivariance", xlim=c(0, 550), ylim=c(0.0, 0.25))

Figure 4: Exponential model fit on Total_C.




 

3 Compute undirected (isotropic) empirical semivariogram for 0.5-m nested points


 

3.1 Compute undirected (isotropic) empirical semivariogram
library(gstat)
## [1] "SpatialPointsDataFrame"
## attr(,"package")
## [1] "sp"
var1 <- variogram(z ~ 1, ~x.axis + y.axis, cutoff = 10, width = 0.5, data = JY)

# create figure of empirical semivariogram
plot(var1, plot.numbers = TRUE, xlab = "Distance", ylab = "Semivariance", cex = 1, cex.axis = 2, xlim=c(0, 3), ylim=c(0.0, 0.25))

Figure 3: Undirected semivariogram of Total_C for points seperated by 0.5-m.






** SUMMARY STATISTICS SKEET SHOOT FIELD**


4. Summary statistics of sampling distance and soil effects on Total_C:

Sampling.dist # Soil series and sampling distance (0.5; 15 and 30 meters) as dependent variables
##                    Df Sum Sq Mean Sq F value   Pr(>F)    
## distance            2  0.286  0.1428   1.200 0.304345    
## S_Series            2  1.969  0.9843   8.271 0.000401 ***
## distance:S_Series   4  1.296  0.3239   2.722 0.031969 *  
## Residuals         141 16.781  0.1190                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1


Table 1:

Descriptive statistics of soil organic C (Total_C) concentration, nitrogen concentration (Total_N), and organic matter (OM). Values are based on the 1st sample set (n = 150). Soil properties represent the top 15 cm of Duke Farm’s grassland (skeet shoot field) soil

N C OM
Soil series n Min Median Max Mean SD Min Median Max Mean SD Min Median Max Mean SD
KkoC 60 0.108 0.170 0.227 0.172a 0.030 0.975 1.684 2.360 1.673a 0.328 1.681 2.903 4.069 2.885a 0.565
PenB 65 0.094 0.161 0.310 0.161a 0.034 0.767 1.556 3.251 1.591a 0.403 1.322 2.683 5.605 2.744a 0.695
PeoC 25 0.099 0.141 0.189 0.144b 0.021 0.853 1.364 1.849 1.347b 0.266 1.471 2.352 3.188 2.323b 0.458

where N = Total Nitrogen, C = Total Carbon and OM = organic matter



Table 2:

Analysis of variance for soil organic C (Total_C) concentration, nitrogen concentration (Total_N), and organic matter (OM) for Duke Farm’s grassland (skeet shoot field) soil n = 150.

Parameters df Sum Sq Mean Sq F value Pr(>F)
Total_C S_Series 2 1.881 0.940 7.492 0.001***
Residuals 147 18.450 0.126 - -
Total_N S_Series 2 0.014 0.007 7.481 0.001***
Residuals 147 0.137 0.001 - -
OM S_Series 2 5.590 2.795 7.492 0.001***
Residuals 147 54.845 0.373 - - 
model2<-aov(Total_C ~ S_Series, data=mcPS)
AoVTotal_C <- summary(model2)

How table 2 was obtained:

AoVTotal_C # Soil series and Total C
##              Df Sum Sq Mean Sq F value   Pr(>F)    
## S_Series      2  1.881  0.9403   7.492 0.000798 ***
## Residuals   147 18.450  0.1255                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1






4.1. Spatial plot of samples:

You can use the plotting functions spplot or bubble as illustrated below (note: the x- and y-axis are the spatial coordinates)**

data = data.frame(x,y,Total_C)
library(sp)
coordinates(data) = ~x+y
class(data)
## [1] "SpatialPointsDataFrame"
## attr(,"package")
## [1] "sp"
spplot(data, "Total_C", colorkey = TRUE, main = "30m Triangular sampling C concentrations (n = 150 Random)")






4.2. IMPACT OF SOIL SERIES ON CARBON CONTENT:

## 
## Study: IMPACT OF SOIL SERIES ON CARBON CONTENT
## 
## HSD Test for Total_C 
## 
## Mean Square Error:  0.1190116 
## 
## S_Series,  means
## 
##       Total_C       std  r   Min   Max
## KkoC 1.673217 0.3278467 60 0.975 2.360
## PenB 1.591369 0.4033851 65 0.767 3.251
## PeoC 1.347400 0.2657072 25 0.853 1.849
## 
## Alpha: 0.05 ; DF Error: 141 
## Critical Value of Studentized Range: 3.349881 
## 
## Groups according to probability of means differences and alpha level( 0.05 )
## 
## Treatments with the same letter are not significantly different.
## 
##       Total_C groups
## KkoC 1.673217      a
## PenB 1.591369      a
## PeoC 1.347400      b

Post-Hoc Test:

plot(out)

4.3. IMPACT OF SAMPLING DISTANCE ON CARBON CONTENT:

## 
## Study: IMPACT OF SAMPLING DISTANCE ON CARBON CONTENT
## 
## HSD Test for Total_C 
## 
## Mean Square Error:  0.1190116 
## 
## distance,  means
## 
##      Total_C       std  r   Min   Max
## 0.5 1.636875 0.4084768 56 0.975 2.764
## 15  1.523259 0.3120135 27 0.874 2.200
## 30  1.563045 0.3554102 67 0.767 3.251
## 
## Alpha: 0.05 ; DF Error: 141 
## Critical Value of Studentized Range: 3.349881 
## 
## Groups according to probability of means differences and alpha level( 0.05 )
## 
## Treatments with the same letter are not significantly different.
## 
##      Total_C groups
## 0.5 1.636875      a
## 30  1.563045      a
## 15  1.523259      a

Post-Hoc Test:

plot(out)




# The End


 

[Back to Top]



Continue to Grazed Fields Data



Continue to Land use effect Data