Topic: Spatial
Topic Description:
Functions for calculating Ripley's K from CTFS R Analytical Tables, many routines for quadrat-based calculations. and calculations of wavelet variance.
File: spatial/block.analysis.r
| View File Source | Download File | No help file available |
Function: wavelet.allsp
Function Description: wavelet.allsp
Function to calculate the wavelet variance curves for all species in one plot. Requires spatstat. Author: Matteo Detto and Tania BrenesOutput is a matrix with the wavelet curve per species called $variances.
Function Arguments:
| Argument | Default Value |
|---|---|
| censdata | |
| plotdim | c(1000,500) |
| gridsize | 2.5 |
| mindbh | NULL |
Arguments Description:
censdata (): census data for the plot containing the variables gx, gy, dbh, status, and sp code;
plotdim c(1000,500): vector with two numbers indicating the plot size;
Sample Usage:
load("bci.full1.rdata")
install.packages('spatastat') library(spatstat) bci.wavelets = wavelet.allsp(bci.full7, plotdim=c(1000,500), gridsize=5)
Function Source:
#
#
# wavelet.allsp
#
#
# Function to calculate the wavelet variance curves for all species in one plot. Requires spatstat.
# Author: Matteo Detto and Tania Brenes
# Output is a matrix with the wavelet curve per species called $variances.
#
#
# censdata (): census data for the plot containing the variables gx, gy, dbh, status, and sp code;
# plotdim c(1000,500): vector with two numbers indicating the plot size;
#
#
# load("bci.full1.rdata")
# install.packages('spatastat')
# library(spatstat)
# bci.wavelets = wavelet.allsp(bci.full7, plotdim=c(1000,500), gridsize=5)
#
#
wavelet.allsp = function(censdata, plotdim=c(1000,500), gridsize=2.5, mindbh=NULL)
{
ptm <- proc.time()
if (is.null(mindbh)) censdata = subset(censdata, status=="A" & !is.na(gx) & !is.na(gy)) else
censdata = subset(censdata, status=="A" & !is.na(gx) & !is.na(gy), dbh>=mindbh)
censdata$sp = factor(censdata$sp)
poly = owin(xrange=c(0,plotdim[1]), yrange=c(0,plotdim[2]))
censdata.ppp = ppp(censdata$gx, censdata$gy, marks=censdata$sp, window=poly)
splitdata = split(censdata.ppp, as.factor(censdata$sp))
counts = quadratcount(splitdata, xbreaks=seq(0, plotdim[1], gridsize), ybreaks= seq(0, plotdim[2], g
ridsize))
n = length(splitdata)
variance = numeric()
sp.density = matrix(NA, ncol=2, nrow=n)
dimnames(sp.density) = list(names(counts), c("number", "density"))
for (i in 1:n)
{
x = wavelet.var(counts[[i]], plotdim=plotdim, gridsize=gridsize, type='raster', k0=8, dj=0.15, gra
ph=FALSE) variance = rbind(variance, x$variance)
sp.density[i,1] = npoints(splitdata[[i]])
sp.density[i,2] = npoints(splitdata[[i]])/(plotdim[1]*plotdim[2])
if (i %in% seq(10,n+10,10)) cat( i, "of", n, " elapsed time = ", (proc.time()-ptm)[3]/60, "minute
s" , "\n") }
dimnames(variance) <- list(names(splitdata), paste("scale",1:ncol(variance)))
cat( "Total elapsed time = ", (proc.time()-ptm)[3]/60, "minutes" , "\n")
return(list(scale = x$scale, variance = variance, density= sp.density, plotdim=plotdim, gridsize=gri
dsize, norm=x$norm, UCL=x$UCL, LCL=x$LCL)) }
#
# wavelet.allsp
#
#
# Function to calculate the wavelet variance curves for all species in one plot. Requires spatstat.
# Author: Matteo Detto and Tania Brenes
# Output is a matrix with the wavelet curve per species called $variances.
#
#
# censdata (): census data for the plot containing the variables gx, gy, dbh, status, and sp code;
# plotdim c(1000,500): vector with two numbers indicating the plot size;
#
#
# load("bci.full1.rdata")
# install.packages('spatastat')
# library(spatstat)
# bci.wavelets = wavelet.allsp(bci.full7, plotdim=c(1000,500), gridsize=5)
#
#
wavelet.allsp = function(censdata, plotdim=c(1000,500), gridsize=2.5, mindbh=NULL)
{
ptm <- proc.time()
if (is.null(mindbh)) censdata = subset(censdata, status=="A" & !is.na(gx) & !is.na(gy)) else
censdata = subset(censdata, status=="A" & !is.na(gx) & !is.na(gy), dbh>=mindbh)
censdata$sp = factor(censdata$sp)
poly = owin(xrange=c(0,plotdim[1]), yrange=c(0,plotdim[2]))
censdata.ppp = ppp(censdata$gx, censdata$gy, marks=censdata$sp, window=poly)
splitdata = split(censdata.ppp, as.factor(censdata$sp))
counts = quadratcount(splitdata, xbreaks=seq(0, plotdim[1], gridsize), ybreaks= seq(0, plotdim[2], g
ridsize))
n = length(splitdata)
variance = numeric()
sp.density = matrix(NA, ncol=2, nrow=n)
dimnames(sp.density) = list(names(counts), c("number", "density"))
for (i in 1:n)
{
x = wavelet.var(counts[[i]], plotdim=plotdim, gridsize=gridsize, type='raster', k0=8, dj=0.15, gra
ph=FALSE) variance = rbind(variance, x$variance)
sp.density[i,1] = npoints(splitdata[[i]])
sp.density[i,2] = npoints(splitdata[[i]])/(plotdim[1]*plotdim[2])
if (i %in% seq(10,n+10,10)) cat( i, "of", n, " elapsed time = ", (proc.time()-ptm)[3]/60, "minute
s" , "\n") }
dimnames(variance) <- list(names(splitdata), paste("scale",1:ncol(variance)))
cat( "Total elapsed time = ", (proc.time()-ptm)[3]/60, "minutes" , "\n")
return(list(scale = x$scale, variance = variance, density= sp.density, plotdim=plotdim, gridsize=gri
dsize, norm=x$norm, UCL=x$UCL, LCL=x$LCL)) }