Skip to content

Commit

Permalink
Update wDist, man
Browse files Browse the repository at this point in the history
  • Loading branch information
bumbanian committed Jan 28, 2022
1 parent e57930f commit daf0dd9
Show file tree
Hide file tree
Showing 3 changed files with 29 additions and 14 deletions.
18 changes: 15 additions & 3 deletions R/wDist.R
Original file line number Diff line number Diff line change
@@ -1,4 +1,4 @@
wDist = function(pdR, sites){
wDist = function(pdR, sites, maxpts = 1e5){

if(class(pdR) != "RasterLayer" & class(pdR) != "RasterStack" & class(pdR) != "RasterBrick"){
stop("input probability density map (pdR) should be one of the following class: RasterLayer, RasterStack or RasterBrick")
Expand All @@ -22,6 +22,13 @@ wDist = function(pdR, sites){
stop("sites should be a SpatialPoints object")
}

if(!is.numeric(maxpts)){
stop("maxpts must be numeric")
}
if(!(round(maxpts) == maxpts) | maxpts < 1){
stop("maxpts must be a positive integer")
}

#make space
wd = list()

Expand All @@ -31,8 +38,13 @@ wDist = function(pdR, sites){

for(i in seq_along(sites)){
pdSP = rasterToPoints(pdR[[i]], spatial = TRUE)
if(length(pdSP) > maxpts){
index = sample(seq(length(pdSP)), maxpts)
pdSP = pdSP[index,]
pdSP@data[,1] = pdSP@data[,1] / sum(pdSP@data[,1])
}

#for safety, using projected data works on most platforms
#for safety; using projected data works on most platforms
pdSP = spTransform(pdSP, "+proj=longlat +ellps=WGS84")
sites = spTransform(sites, "+proj=longlat +ellps=WGS84")

Expand Down Expand Up @@ -280,7 +292,7 @@ plot.wDist = function(x, ..., bin = 20, pty = "both", index = c(1:5)){
}
for(j in seq_along(bins)){
xy = wedge(bins[j], bins[j] + bin, vals[j])
c = col2rgb(index[i])
c = col2rgb(i)
polygon(xy, col = rgb(c[1], c[2], c[3],
alpha = 200, maxColorValue = 255))
}
Expand Down
4 changes: 2 additions & 2 deletions man/getIsoscapes.Rd
Original file line number Diff line number Diff line change
Expand Up @@ -29,10 +29,10 @@ Accepted \code{isoType} values are:
\item{"GlobalPrecipGS"}{Global growing-season precipitation H and O isotope values}
\item{"GlobalPrecipMA"}{Global mean-annual precipitation H and O isotope values}
\item{"GlobalPrecipMO"}{Global monthly precipitation H and O isotope values}
\item{"GlobalPrecipMA"}{Global mean-annual and monthly precipitation H and O isotope values}
\item{"GlobalPrecipALL"}{Global mean-annual and monthly precipitation H and O isotope values}
\item{"USPrecipMA"}{High-resolution contiguous USA mean-annual precipitation H and O isotope values}
\item{"USPrecipMO"}{High-resolution contiguous USA monthly precipitation H and O isotope values}
\item{"USPrecipMA"}{High-resolution contiguous USA mean-annual and monthly precipitation H and O isotope values}
\item{"USPrecipALL"}{High-resolution contiguous USA mean-annual and monthly precipitation H and O isotope values}
\item{"USSurf"}{High-resolution contiguous USA surface water H and O isotope values}
\item{"USTap"}{High-resolution contiguous USA surface water H and O isotope values}
\item{"USSr"}{High-resolution contiguous USA Sr isotope ratios}
Expand Down
21 changes: 12 additions & 9 deletions man/wDist.Rd
Original file line number Diff line number Diff line change
Expand Up @@ -11,22 +11,25 @@ Calculate the distance and bearing of migration for one or more samples, weighte
}

\usage{
wDist(pdR, sites)
wDist(pdR, sites, maxpts = 1e5)
}

\arguments{
\item{pdR}{
RasterStack or RasterBrick of n probability density maps, e.g., as produced by \code{pdRaster}.
}
RasterStack or RasterBrick of n probability density maps, e.g., as produced by \code{pdRaster}.
}
\item{sites}{
SpatialPoints* object containing the collection locations for the n samples represented in pdR.
}
SpatialPoints* object containing the collection locations for the n samples represented in \code{pdR}.
}
\item{maxpts}{
numeric. Maximum number of grid cells at which to calculate bearing and distance.
}
}

\details{
\code{pdR} and \code{sites} must be of equal length and corresponding order. Names in the output object are taken from the names of the layers in \code{pdR}.

Distances and bearings are calculated on the WGS84 geoid using functions from the \pkg{geosphere} package.
Distances and bearings are calculated on the WGS84 geoid using functions from the \pkg{geosphere} package. These calculations can take a long time for large rasters. If \code{maxpts} is less than the number of grid cells in each \code{pdR} layer, calculations are caried out for \code{maxpts} randomly selected cells.

Bearing values correspond to the initial bearing from source to collection location, and are reported on a scale of -180 to +180 degrees. The statistical metrics are rectified so that values for distributions spanning due south are reported correctly. Both weighted bearing and distance distributions are often multimodal, and it is recommended to review the distribution densities to assess the representativeness of the statistics (e.g., using \code{\link{plot.wDist}}).
}
Expand All @@ -35,10 +38,10 @@ Bearing values correspond to the initial bearing from source to collection locat
Returns an object of class \dQuote{wDist}, a list of length n. Each item contains three named objects:
\item{stats}{
named number. Statistics for the distance (dist, meters) and bearing (bear, degrees) between source and collection locations, including the weighted mean (wMean) and quantile (wXX) values.}
\item{d.density}{
\item{d.dens}{
density. Weighted kernel density for the distance between source and collection locations (meters). See \code{\link[stats]{density}}.
}
\item{s.density}{
\item{b.dens}{
density. Weighted kernel density for the bearing between source and collection locations (degrees). See \code{\link[stats]{density}}.
}
}
Expand Down Expand Up @@ -71,6 +74,6 @@ sites = d$data[sample(seq(length(d$data)), 4),]
wd = wDist(pp, sites)

# structure of the wDist object
str(wd[[1]])
str(wd, 2)
}

0 comments on commit daf0dd9

Please sign in to comment.