In NINAnor/SurveyPower: Functions for calculating power of survey design
options(list("xtable.include.rownames" = F,
"xtable.comment" = F))
require(knitr)
knitr::opts_chunk$set(echo = TRUE, tidy = 'styler',
dev = c("png", "pdf", "svg"),
dev.args=list(family='Arial'),
dpi = 600, fig.path = "figure/")
require(extrafont)
#loadfonts()
Intro
We here explore some effects of rare occurrences and imperfect detection in observation of species, and explore alternative survey regimes. This is mostly relevant for rare species, and not for the general survey of insects. A use case is early detection of alien species, which we explore a little deeper.
The objective here is to identify suitable sampling strategies for a couple of different scenarios, where we maximise the chances of detecting a species at least once while considering the economic costs. The basic principle here is that a species can be present in a sampling location with a given probability of occurrence. The observers then have a given chance (probability) of detecting a species on a visit, if it is present. The resulting probability of detecting a species will thus be dependent on both the probability of occurrence and the probability of detection, how many locations are visited, and how many visits we make at each location.
If $\psi$ represents the occurrence probability and $\theta$ is the detection probability, the probability of observing the species at least once in $J$ locations visited $K$ times is: $$probObserv = 1 - (1 - \psi * (1 - (1 - \theta)^K))^J.$$ To simplify, we here assume that the costs scale roughly linearly to the total amount of visits, i.e. $totalCost = J * K * c$, where $c$ is,the cost of one visit/survey.
First look
We can explore how the overall probability of observing the species, and the associated costs depend on $\psi, \theta, J, K, c$ graphically. I have made a convenience function obsProb to calculate the values. A lot of the results seem pretty obvious in retrospect, but it is still worth plotting them to get a better feel for the possibilities. At the moment, we only consider one and the same cost for each survey visit. It may be worth while coding up using a higher initial cost of the first visit, and lower costs for subsequent visits, to account for the economics of scale.
require(SurveyPower)
##custom library, available through devtools::install_github("NINAnor/SurveyPower")
require(dplyr)
require(gridExtra)
require(ggplot2)
require(xtable)
We can specify a single or a range of values for the occurrence probability in each location, the detection probabilities, the number of locations visited, the number of visits per location, and the cost of each visit. To illustrate, we here use occurrence probabilities ranging from 0.01 to 1, two different detection probabilities as 0.5, and 0.8. We specify 30 locations, each visited 4 times, at an individual visit cost of 5000.
obsDf <- obsProb(occProb = seq(0.01, 1, by = 0.05),
detectProb = c(0.5, 0.8),
locations = 30,
visits = 1,
visitCost = 5000)
obsDf
The probability of observing the species at least once is found in the obsProb column, and the total cost in the column totCost. The function returns an object of a specific class with a custom plotting function (see figure \ref{princPlot}). This makes repeated plottings easier. We can specify a grouping variable to split up the lines.
```r"}
plot(obsDf, group = "detectProb",
xVar = "occProb",
yVar = "obsProb",
titleVar = c("locations", "visits"),
hline = 0.8)
Figure \ref{princPlot} shows how the overall obervation probability is dependent on both occurrence and detection probability. A threshold of a total observation probability of 0.8 is added for comparison. The threshold is reached in this case when the probability of occurrence in each location is 6%. With 30 locations, the overall observation probability rises quite sharply as a result of increased occurrence probability. With only 4 visits per location, the detection probability limits the overall achievable observation probability.
However, these occurrence probabilities are probably unreasonably high for rare species in the real world. If we assume that we search for a truly rare species with an occurrence probability of only 0.001, we find that reaching overall detectabilities above 80% is challenging (Figure \ref{rareLocPlot}).
```r
rareLocDf <- obsProb(occProb = 0.001,
detectProb = seq(0.4, 0.8, by = 0.2),
locations = seq(30, 1000, by = 20),
visits = seq(4, 16, by = 4),
visitCost = 5000)
rareLocDf
```r"}
plot(rareLocDf, group = "visits",
xVar = "locations",
yVar = "obsProb",
hline = 0.8) +
facet_grid(occProb ~ detectProb)
In figure \ref{rareLocPlot}, we see that observing a very rare species with a high certainty is difficult even with a very large number of visited locations. In these cases, it doesn't really help to visit each location many times, as the overall probability is limited by the number of locations we visit. Figure \ref{rareVisitPlot} shows the results of maximising the number of visits in fixed, but large number of locations, for a very rare species.
```r
rareVisitDf <- obsProb(occProb = 0.001,
detectProb = seq(0.2, 0.8, by = .2),
locations = 250,
visits = seq(1, 11, by = 5),
visitCost = 5000)
rareVisitDf
```r"}
plot(rareVisitDf, group = "detectProb",
xVar = "visits",
yVar = "obsProb",
titleVar = c("occProb", "locations"),
hline = 0.8)
Although the overall cost increases linearly with the total number of samples (figure \ref{rareCostPlot}), in cases with very rare species, this doesn't mean that the overall detection probability continues to increase indefinitely (figure \ref{rareVisitPlot})
```r
rareCostDf <- obsProb(occProb = 0.001,
detectProb = 0.4,
locations = seq(50, 250, by = 50),
visits = seq(1, 21, by = 5),
visitCost = 5000)
rareCostDf
```r"}
plot(rareCostDf, group = "locations",
xVar = "visits",
yVar = "totCost",
titleVar = "visitCost")
<!-- ```r"} -->
<!-- plot(rareVisitDf, group = "detectProb", -->
<!-- xVar = "totCost", -->
<!-- yVar = "obsProb", -->
<!-- titleVar = c("occProb", "locations")) -->
<!-- ``` -->
Sampling strategies
==============
When we deal with a very rare species, it is little use increasing the number of visits to each location (or to spend money maximizing the detection probability), if we only visit a few locations. In these cases, we have to cover a very large number of locations, just to have a chance of visiting a location with a presence (figure \ref{lowOccurrPlot}). Still, for a very rare species such as displayed in figure \ref{lowOccurrPlot}, with an occurrence probability of 0.01%, reaching an overall observation probability of 80% requires more at least 1650 locations, which would be unfeasible for most survey programs.
Alternatively, if detectability is low but presence is relatively high, we should focus on increasing the number of revisits per location, instead of trying to cover many locations (figure \ref{lowDetectPlot}).
```r
lowOccurrDf <- obsProb(occProb = 0.001,
detectProb = 0.8,
locations = seq(50, 2000, by = 50),
visits = c(1, seq(5, 20, by = 5)),
visitCost = 5000)
lowOccurrDf
```r"}
plot(lowOccurrDf, group = "visits",
xVar = "locations",
yVar = "obsProb",
titleVar = c("occProb", "detectProb"),
hline = 0.8)
```r
lowDetectDf <- obsProb(occProb = 0.05,
detectProb = 0.2,
locations = seq(50, 250, by = 50),
visits = c(1, seq(5, 20, by = 5)),
visitCost = 5000)
lowDetectDf
```r"}
plot(lowDetectDf, group = "visits",
xVar = "locations",
yVar = "obsProb",
titleVar = c("occProb", "detectProb"),
hline = 0.8)
Plausible values for early detection of alien species.
-------------
If we deal with for example early detection of alien species, it is reasonable to assume that we only have to consider rather low occurrence probabilities. We can explore our possibilities of observing a species that occurr in between 0.1 to 1% of all studied locations. For detection probabilities, we choose a range of expected, while hopeful, values. We her set the number of locations to 100 and with two visits, as a reasonable setup to complete within a year.
```r
guessDf <- obsProb(occProb = c(0.001, 0.005, 0.01, 0.02, 0.04, 0.05, 0.1),
detectProb = c(0.1, seq(0.2, 0.8, by = 0.2)),
locations = c(50, seq(100, 300, by = 100)),
visits = seq(2, 6, by = 2),
visitCost = 5000)
guessDf
plot(guessDf, group = "detectProb",
xVar = "occProb",
yVar = "obsProb",
hline = 0.8) +
facet_grid(visits ~ locations) +
theme(axis.text.x = element_text(size = 5))
```r"}
plot(guessDf, group = "detectProb",
xVar = "occProb",
yVar = "obsProb",
hline = 0.8) +
facet_grid(visits ~ locations) +
scale_y_log10() +
scale_x_log10() +
theme(axis.text.x = element_text(size = 6)) +
ylab("Observasjonssannsynlighet") +
xlab("Forekomstsannsynlighet") +
guides(color=guide_legend(title="Pdeteksjon"))
Unequally distributed probabilities
==================
So far, we have only considered equally distributed probabilities for occurrences, meaning that every grid cell has equal probability of containing an occurrence. In the equation $probObserv = (1 - (1 - \psi)^J) * (1 - (1 - \theta)^K))$, $\psi$ designates the probability that a site we visit contains a certain species. This is the parameter which up until now have been the same for all sites. This is not such a good model for many real world situations. Often we have some knowledge, or at least guesses of where a species is more likely to be present. Therefore, the probability that a specific site that you visit contains a certain species will vary between sites. It will depend both on the probability of a site containing a species, and the probability that you visit that site.
We can view the probability of a specific site containing a species as a weighted sampling without replacement. For example, if we know there are a 100 sites containing species x, we can calculate the probability that each site contains species x if we know the probability weights. This probability can be written, following @Efraimidis_2006, as $p_i(k) = \frac{w_i} {\sum_{S_j \in V - S} {w_j}}$. For brevity, however, we will simply designate this probability as $Pr[w_i, n]$, where the weights $w_i$ sums to 1, and $n$ is the number of sites with the species. Similarly, we choose which sites to visit as a sampling without replacement, so that the probability of visiting a site $i$ is $Pr[u_i, v]$, where $u_i$ is the visitation weights, which sums to 1, and $v$ is the number of sites you visit. The probability of visiting a site with an alien species then becomes $\psi_i = Pr[w_i, n] * Pr[u_i, v]$, and the probability of visiting any location inhabited by an alien species is $PoccurrVisit = (1 - \prod_{i}^{}(1 - (Pr[w_i, n] * Pr[u_i, v])))$.
For simplicity, we can assume that the probability of detection is equal for each site and visit. The probability of detecting an alien species is then $Pdetect = (1 - \prod_{i}^{}(1 - (Pr[w_i, n] * Pr[u_i, v]))) * (1 - (1 - d)^K$, where d is the detection probability, and K is again the number of visits per site.
In the best of worlds, $w_i$ and $u_i$ would match up well, so that we visit the most likely sites to contain a species of interest. The overall probability of detection will decrease as the difference between these variables increase, or in simpler terms we visit the wrong places. Also, the probability will go down the more spread out the probability weights are. In the extreme case, with no spread in these probabilities, there is 100 % certainty that a species will be present in location 1, and 100 % certainty that we will visit just that. In this case, the probabilities are 1 for both occurrence and visitation, so that we are certain we visit a site with a presence. In the other end, there might be no information in the weights, so that the occurrences are randomly spread out, and we visit random sites. In that case, $w_i$ are all the same and $u_i$ are all the same, and we end up with the first equation.
In reality, we don't know the true weights that the sites will contain a specific species, and we might choose to visit the wrong sites. For this example calculation, we will assume we know the occurrence weights, and visit them accordingly. In other words, that $w_i = u_i$. If we stipulate the desired detection probability (at e.g 0.8), we can calculate how many sites $v$ we need to visit to be able to detect a species that occurr in $n$ sites with weighted probabilities $w_i$, with a specific certainty.
So far, the best estimate for the occurrence weights are the occurrence modelling of alien vascular plants from Olsen et al. 2017. Using this as input, and some simple assumptions, we can calculate the needed number of sites.
To get test it out, we can use the 10km scale, which isn't so resource intensive. The 1km scale isn't really interactive since it takes to long time to calculate.
```r
require(RPostgreSQL)
require(DBI)
require(rpostgis)
require(SurveyPower)
require(tidyverse)
data(map10km)
con <- dbConnect(RPostgreSQL::PostgreSQL(), host = "gisdata-db.nina.no", dbname = "gisdata", user = "postgjest", password = "gjestpost")
#pred <- pgGetRast(con, name = c("hotspot_ias", "bigPred1km"), rast = "rast", bands = 1, boundary = NULL)
predQ <- "SELECT ssbid, COALESCE(avg(ST_Value(pred.rast, ST_Centroid(ssb.geom))), 0) pred
FROM ssb_data_utm33n.ssb_10km ssb,
hotspot_ias.\"evenintbigpred1km\" pred
WHERE ST_Intersects(ssb.geom, pred.rast)
GROUP BY ssbid"
pred <- dbGetQuery(con, predQ)
pred <- pred %>%
mutate(ssbid = as.character(ssbid))
predMap <- map10km %>% left_join(pred, by = c("ssbid" = "ssbid"))
predMap$pred[is.na(predMap$pred)] <- 0
plot(predMap["pred"],
border = 0.2)
predWeights10km <- predMap %>%
select(sites = ssbid,
weights = pred) %>%
sf::st_set_geometry(NULL)
devtools::use_data(predWeights10km)
data(predWeights10km)
But we start with a situation where the occupied patches are distributed randomly, and we visit the sites randomly. In other words, where the occurrence and visitation weights all are equal.
predWeightsZero<- predWeights10km
predWeightsZero$weights <- 1
system.time(predProb50Zero <- weightedDetection(occWeights = predWeightsZero,
visWeights = predWeights10km,
noOccur = 5,
noLocations = seq(50, 200, by = 50),
noVisits = 1,
detectProb = 1))
predProb50Zero
plot(predProb50Zero, threshold = 0.8)
We can see that the probability of visiting a randomly occupied cell in this case starts from about 0.07 and approaches 0.2 as we increase our number of visited locations from 50 to 200. We can quality check this with a simpler function.
test <- function(noOccur = 50,
noLocations = 5,
nIter = 999){
prop <- function(noLocations. = noLocations,
noOccur. = noOccur){
visited <- sample(1:4057, noLocations., replace = F) #number of 10km cells
occupied <- sample(1:4057, noOccur, replace = F)
any(visited %in% occupied)
}
sum(replicate(nIter, prop())/nIter)
}
test()
We get the same results if the occurrence probabilities is distributed according to informative weights, and only the visitations are random (not shown).
But what happens when we have information about the occurence of the species? In effect, we limit the number of potential sites we visit to a smaller value, which have higher probability of housing the species. We use the prediction map to set the occurrences, and visitation probabilities. We continue with the detection probability set to 1, with just 1 visit per site.
system.time(predProb <- weightedDetection(occWeights = predWeights10km,
visWeights = predWeights10km,
noOccur = 5,
noLocations = seq(50, 200, by = 50),
noVisits = 1,
detectProb = 1,
nIter = 999))
predProb
plot(predProb, threshold = 0.8)
We then reach quite high probabilities for observing the species. This of course depend on the quality of the predictions, i.e. our weights. This prediction map is actually quite informative. From the histogram of weights, we can see that a small number of sites have a proportionally high weight. This limits the realised occurrences quite a bit.
hist(predWeights10km$weights)
Calculations based on 1km map
require(RPostgreSQL)
require(DBI)
require(rpostgis)
require(SurveyPower)
require(tidyverse)
data(map1km)
con <- dbConnect(RPostgreSQL::PostgreSQL(), host = "gisdata-db.nina.no", dbname = "gisdata", user = "postgjest", password = "gjestpost")
#pred <- pgGetRast(con, name = c("hotspot_ias", "bigPred1km"), rast = "rast", bands = 1, boundary = NULL)
predQ <- "SELECT ssbid, ST_Value(pred.rast, ST_Centroid(ssb.geom)) pred
FROM ssb_data_utm33n.ssb_1km ssb,
hotspot_ias.\"evenintbigpred1km\" pred
WHERE ST_Intersects(ssb.geom, pred.rast)
"
pred <- dbGetQuery(con, predQ)
pred <- pred %>%
mutate(ssbid = as.character(ssbid))
predMap1km <- map1km %>% left_join(pred, by = c("ssbid" = "ssbid"))
predMap1km$pred[is.na(predMap1km$pred)] <- 0
#plot(predMap1km["pred"]) #Slow
predWeights1km <- predMap1km %>%
select(sites = ssbid,
weights = pred) %>%
sf::st_set_geometry(NULL)
devtools::use_data(predWeights1km)
Since we now have about 100 times as many potential sites, we would need to increase the occurrences accordingly. But while 50 out of 4057 10x10km squares sounds reasonable for an "early" detection, multiplying this with 100 yields 5000 locations in a 1x1km grid. This sounds like a lot for an early detection. But it puts the 50 occurrences above into perspective. Surveying 10x10km cells with good detection probability is a tall order.
For the 1x1km analysis, the calculations takes to much time to do on the fly, so we pre-calculate them and load the results.
data("predWeights1km")
It is reasonable to assume that we won't reach a higher detection probability than 0.8 for a single visit, and even that is probably high for a truly novel species. But we can explore the range of occupied sites that are reasonable to handle with such a good detection probability.
data(predOccs)
predOcc500Det0.05 <- weightedDetection(occWeights = predWeights1km,
visWeights = predWeights1km,
noOccur = 500,
noLocations = seq(200, 1000, by = 200),
noVisits = 1,
detectProb = 0.8/16,
nIter = 999)
predOcc100Det0.8 <- weightedDetection(occWeights = predWeights1km,
visWeights = predWeights1km,
noOccur = 100,
noLocations = seq(50, 600, by = 50),
noVisits = 1,
detectProb = 0.8,
nIter = 999)
predOcc500Det0.8 <- weightedDetection(occWeights = predWeights1km,
visWeights = predWeights1km,
noOccur = 500,
noLocations = seq(50, 300, by = 50),
noVisits = 1,
detectProb = 0.8,
nIter = 999)
predOcc200Det0.8 <- weightedDetection(occWeights = predWeights1km,
visWeights = predWeights1km,
noOccur = 200,
noLocations = seq(50, 600, by = 50),
noVisits = 1,
detectProb = 0.8,
nIter = 999)
predOcc100Det0.05Vis1 <- weightedDetection(occWeights = predWeights1km,
visWeights = predWeights1km,
noOccur = 100,
noLocations = seq(200, 1000, by = 200),
noVisits = 1,
detectProb = 0.8/16,
nIter = 999)
predOcc500Det0.5 <- weightedDetection(occWeights = predWeights1km,
visWeights = predWeights1km,
noOccur = 500,
noLocations = seq(200, 1000, by = 200),
noVisits = 1,
detectProb = 0.8/16*10,
nIter = 999)
predOcc500Det0.2 <- weightedDetection(occWeights = predWeights1km,
visWeights = predWeights1km,
noOccur = 500,
noLocations = seq(200, 1000, by = 200),
noVisits = 1,
detectProb = 0.8/16*4,
nIter = 999)
save(predOcc500Det0.05, predOcc100Det0.8, predOcc500Det0.8, predOcc200Det0.8, predOcc100Det0.05Vis1, predOcc500Det0.5, predOcc500Det0.2, file = "pred_occs.Rdata")
#collect a bunch of small files into one
predFiles <- list.files(pattern = "pred*")
lapply(predFiles, load, .GlobalEnv)
predObjectsToStore <-grep("predOcc*", ls(), value = T)
save(list = predObjectsToStore, file = "predOccs.Rdata")
```r "}
plot(predOcc500Det0.8,
threshold = 0.8 ,
xlab = "Antall lokaliteter",
ylab = "Observasjonssannsynlighet")
So for the case of 500 occupied cells, we reach the target observation probability after visiting about 200 sites. In similar fashion, we can explore the case with 200 and 100 occupied cells, respectively.
```r
xtable(predOcc200Det0.8)
```r"}
plot(predOcc200Det0.8,
threshold = 0.8,
xlab = "Antall lokaliteter",
ylab = "Observasjonssannsynlighet")
```r
xtable(predOcc100Det0.8)
```r"}
plot(predOcc100Det0.8,
threshold = 0.8)
In the case of 200 occupied cells, we reach an overall observation probability of 0.8 after about 500 visited sites with an observation probability of 0.8. With only 100 occupied cells, we don't reach the target of 0.8 even after 600 visited locations.
We can see the effect of a lower detection probability.
Surveying 250x250m squares
-------------
In practice, it can be challenging to survey even a 1x1km grid with any respectable observation probability. We might therefore subdivide the squares in smaller units, with the result that our detection probability decreases from simply not covering the place where the species occupies. If a species just occurs in one out of 16 250x250 squares within a 1x1km square and we visit only one such smaller square, our detection probability drops to 1/16 of the former level. How much the detection probability drops of course depend on the aggregation pattern of the species, i.e. how much of the 1x1km cell it is present in. On a tangent, the number of occurrences we consider to be acceptable within an "early detection" framework depends on the scale of the grid cells considered, and the way that the species are aggregated. 500 out of ca 50 000 1x1km cells constitutes occurrences in 1% of the cells. If we accept a 1% occurrence in the about 800 000 250x250 cells, that amounts to 8000 occurrences. As long as these are not extremy aggregated, this could hardly be seen as an early establishment phase.
Case of 500 1x1km cells occupied, but we survey only 250x250m subsquares
---------------
We start with the case of 500 occurrences spread out in the 50 000 1x1km cells, but when we visit just a 16th of these cells once.
```r
plot(predOcc500Det0.05,
threshold = 0.8,
xlab = "Antall lokaliteter",
ylab = "Observasjonssannsynlighet")
It is clear that these conditions does not let us reach the desired detection probability of 0.8 even with a great number of sampled locations.
We can continue to explore the possibilities with a lower total occurrence, for example when a species is present in 100 out of the 50 000 1x1km grid cells.
plot(predOcc100Det0.05Vis1,
threshold = 0.8,
xlab = "Antall lokaliteter",
ylab = "Observasjonssannsynlighet")
This situation leaves us with very slim chances of detecting the species.
We can mitigate these low numbers by surveying the same 1x1km square multiple times. Note that I here assume that the species is stationary in only one place, and we return to different subplots within the 1x1km square, so that we cover an increasing portion of the sample location. The detectability thus increases linearly. For example, if we have a detection probability of 0.8 in a given location, but we only visit 10 1/16th of that location, we get a detection probability of 0.8/1610 = 0.5. If we would instead visit 2 subsquares of 1/16, we would get 0.8/164 = 0.2
plot(predOcc500Det0.5,
threshold = 0.8,
xlab = "Antall lokaliteter",
ylab = "Observasjonssannsynlighet")
plot(predOcc500Det0.2,
threshold = 0.8,
xlab = "Antall lokaliteter",
ylab = "Observasjonssannsynlighet")
Surveys limited to a specific area
One strategy is to limit the surveyed areas to focus the efforts. We may thereby increase our chances to visit a location that is occupied by a rare species, but we also cannot find species that are not present in the area of interest.
This situation can be modelled using the same function weightedDetection, by setting the visitation probability to zero for the locations outside our focus area. Based on the prediction map, I have selected a set of kommunes around Oslo that seems to encompass the areas with highest predictions.
data(surveyCalc)
Oslo region
occWeights <- predWeights1km
osloArea <- c("Vestby",
"Svelvik",
"Lørenskog",
"Frogn",
"Ås",
"Skedsmo",
"Oppegård",
"Ski",
"Hobøl",
"Spydeberg",
"Sande",
"Oslo",
"Bærum",
"Asker",
"Lier",
"Øvre Eiker",
"Nedre Eiker",
"Drammen",
"Hof",
"Holmestrand",
"Nesodden",
"Frogn",
"Røyken"
)
osloAreaKommuneNr <- c("0211",
"0711",
"0230",
"0215",
"0214",
"0231",
"0217",
"0213",
"0138",
"0123",
"0713",
"0301",
"0219",
"0220",
"0626",
"0624",
"0625",
"0602",
"0714",
"0702",
"0216",
"0215",
"0627"
)
ssbidToVisit <- map1km %>%
filter(KOMMUNENUMMER %in% osloAreaKommuneNr) %>%
select(ssbid, kommune)
plot(ssbidToVisit["kommune"], key.width = lcm(4), key.pos = 4,
border = 0)
osloVisWeights <- occWeights
osloVisWeights$weights[!(osloVisWeights$sites %in% ssbidToVisit$ssbid)] <- 0
summary(osloVisWeights$weights)
summary(occWeights$weights)
Bergen
One of several possible demarcations of "Bergen area"
bergenKommuneNr <- c("1265",
"1259",
"1253",
"1260",
"1263",
"1251",
"1411",
"1252",
"1256",
"1266",
"1264",
"1242",
"1247",
"1201",
"1243",
"1246",
"1245"
)
ssbidBergen <- map1km %>%
filter(KOMMUNENUMMER %in% bergenKommuneNr) %>%
select(ssbid, kommune)
plot(ssbidBergen["kommune"], key.width = lcm(4), key.pos = 4,
border = 0)
bergenVisWeights <- occWeights
bergenVisWeights$weights[!(bergenVisWeights$sites %in% ssbidBergen$ssbid)] <- 0
summary(bergenVisWeights$weights)
summary(occWeights$weights)
Trondheim
trondheimKommuneNr <- c("1638",
"1657",
"1636",
"1664",
"1662",
"1601",
"1663",
"1714",
"1653")
ssbidTrondheim <- map1km %>%
filter(KOMMUNENUMMER %in% trondheimKommuneNr) %>%
select(ssbid, kommune)
plot(ssbidTrondheim["kommune"],
key.width = lcm(4), key.pos = 4,
border = 0)
trondheimVisWeights <- occWeights
trondheimVisWeights$weights[!(trondheimVisWeights$sites %in% ssbidTrondheim$ssbid)] <- 0
summary(trondheimVisWeights$weights)
summary(occWeights$weights)
Stavanger
occWeights <- predWeights1km
stavangerArea <- c("Hå",
"Time",
"Sandnes",
"Randaberg",
"Stavanger",
"Sola",
"Klepp"
)
stavangerAreaKommuneNr <- c("1119",
"1121",
"1102",
"1127",
"1103",
"1124",
"1120"
)
ssbidToVisit <- map1km %>%
filter(KOMMUNENUMMER %in% stavangerAreaKommuneNr) %>%
select(ssbid, kommune)
plot(ssbidToVisit["kommune"], key.width = lcm(4), key.pos = 4,
border = 0)
stavangerVisWeights <- occWeights
stavangerVisWeights$weights[!(stavangerVisWeights$sites %in% ssbidToVisit$ssbid)] <- 0
summary(stavangerVisWeights$weights)
sum(stavangerVisWeights$weights >0)
Run calculations
oslo47Sites <- weightedDetection(occWeights = occWeights,
visWeights = osloVisWeights,
noOccur = 500,
noLocations = 47,
detectProb = 0.426,
nIter = 999,
noVisits = 1)
oslo47Sites16 <- weightedDetection(occWeights = occWeights,
visWeights = osloVisWeights,
noOccur = 500,
noLocations = 47,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
oslo100Sites16 <- weightedDetection(occWeights = occWeights,
visWeights = osloVisWeights,
noOccur = 500,
noLocations = 100,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
oslo66Sites16 <- weightedDetection(occWeights = occWeights,
visWeights = osloVisWeights,
noOccur = 500,
noLocations = 66,
detectProb = 0.426*1.5/16,
nIter = 999,
noVisits = 1)
oslo200Sites16 <- weightedDetection(occWeights = occWeights,
visWeights = osloVisWeights,
noOccur = 500,
noLocations = 200,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
oslo140Sites16 <- weightedDetection(occWeights = occWeights,
visWeights = osloVisWeights,
noOccur = 500,
noLocations = 140,
detectProb = 0.426*1.5/16,
nIter = 999,
noVisits = 1)
##Bergen
bergen47Sites <- weightedDetection(occWeights = occWeights,
visWeights = bergenVisWeights,
noOccur = 500,
noLocations = 47,
detectProb = 0.426,
nIter = 999,
noVisits = 1)
bergen47Sites16 <- weightedDetection(occWeights = occWeights,
visWeights = bergenVisWeights,
noOccur = 500,
noLocations = 47,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
bergen100Sites16 <- weightedDetection(occWeights = occWeights,
visWeights = bergenVisWeights,
noOccur = 500,
noLocations = 100,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
bergen66Sites16 <- weightedDetection(occWeights = occWeights,
visWeights = bergenVisWeights,
noOccur = 500,
noLocations = 66,
detectProb = 0.426*1.5/16,
nIter = 999,
noVisits = 1)
bergen200Sites16 <- weightedDetection(occWeights = occWeights,
visWeights = bergenVisWeights,
noOccur = 500,
noLocations = 200,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
bergen140Sites16 <- weightedDetection(occWeights = occWeights,
visWeights = bergenVisWeights,
noOccur = 500,
noLocations = 140,
detectProb = 0.426*1.5/16,
nIter = 999,
noVisits = 1)
##Trondheim
trondheim47Sites <- weightedDetection(occWeights = occWeights,
visWeights = trondheimVisWeights,
noOccur = 500,
noLocations = 47,
detectProb = 0.426,
nIter = 999,
noVisits = 1)
trondheim47Sites16 <- weightedDetection(occWeights = occWeights,
visWeights = trondheimVisWeights,
noOccur = 500,
noLocations = 47,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
trondheim100Sites16 <- weightedDetection(occWeights = occWeights,
visWeights = trondheimVisWeights,
noOccur = 500,
noLocations = 100,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
trondheim66Sites16 <- weightedDetection(occWeights = occWeights,
visWeights = trondheimVisWeights,
noOccur = 500,
noLocations = 66,
detectProb = 0.426*1.5/16,
nIter = 999,
noVisits = 1)
trondheim200Sites16 <- weightedDetection(occWeights = occWeights,
visWeights = trondheimVisWeights,
noOccur = 500,
noLocations = 200,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
trondheim140Sites16 <- weightedDetection(occWeights = occWeights,
visWeights = trondheimVisWeights,
noOccur = 500,
noLocations = 140,
detectProb = 0.426*1.5/16,
nIter = 999,
noVisits = 1)
norge140Sites16 <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = 140,
detectProb = 0.426*1.5/16,
nIter = 999,
noVisits = 1)
#Norge
norge47Sites <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = 47,
detectProb = 0.426,
nIter = 999,
noVisits = 1)
norge47Sites16 <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = 47,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
norge100Sites16 <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = 100,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
norge66Sites16 <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = 66,
detectProb = 0.426*1.5/16,
nIter = 999,
noVisits = 1)
norge200Sites16 <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = 200,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
##small expensive
#2 visits, 1.5 detection prob
norge50Sites16Vis2 <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = 50,
detectProb = 0.426*1.5/16*2,
nIter = 999,
noVisits = 1)
oslo50Sites16Vis2 <- weightedDetection(occWeights = occWeights,
visWeights = osloVisWeights,
noOccur = 500,
noLocations = 50,
detectProb = 0.426*1.5/16*2,
nIter = 999,
noVisits = 1)
bergen50Sites16Vis2 <- weightedDetection(occWeights = occWeights,
visWeights = bergenVisWeights,
noOccur = 500,
noLocations = 50,
detectProb = 0.426*1.5/16*2,
nIter = 999,
noVisits = 1)
trondheim50Sites16Vis2 <- weightedDetection(occWeights = occWeights,
visWeights = trondheimVisWeights,
noOccur = 500,
noLocations = 50,
detectProb = 0.426*1.5/16*2,
nIter = 999,
noVisits = 1)
save(oslo47Sites, oslo47Sites16, oslo66Sites16, oslo100Sites16, oslo140Sites16, oslo200Sites16, oslo50Sites16Vis2,
bergen47Sites, bergen47Sites16, bergen66Sites16, bergen100Sites16, bergen140Sites16, bergen200Sites16, bergen50Sites16Vis2,
trondheim47Sites, trondheim47Sites16, trondheim66Sites16, trondheim100Sites16, trondheim140Sites16, trondheim200Sites16, trondheim50Sites16Vis2,
norge47Sites, norge47Sites16, norge66Sites16, norge100Sites16, norge140Sites16, norge200Sites16, norge50Sites16Vis2,
file = "survey_calc.Rdata")
stavanger47Sites <- weightedDetection(occWeights = occWeights,
visWeights = stavangerVisWeights,
noOccur = 500,
noLocations = 47,
detectProb = 0.426,
nIter = 999,
noVisits = 1)
stavanger47Sites16 <- weightedDetection(occWeights = occWeights,
visWeights = stavangerVisWeights,
noOccur = 500,
noLocations = 47,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
stavanger100Sites16 <- weightedDetection(occWeights = occWeights,
visWeights = stavangerVisWeights,
noOccur = 500,
noLocations = 100,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
stavanger66Sites16 <- weightedDetection(occWeights = occWeights,
visWeights = stavangerVisWeights,
noOccur = 500,
noLocations = 66,
detectProb = 0.426*1.5/16,
nIter = 999,
noVisits = 1)
stavanger200Sites16 <- weightedDetection(occWeights = occWeights,
visWeights = stavangerVisWeights,
noOccur = 500,
noLocations = 200,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
stavanger140Sites16 <- weightedDetection(occWeights = occWeights,
visWeights = stavangerVisWeights,
noOccur = 500,
noLocations = 140,
detectProb = 0.426*1.5/16,
nIter = 999,
noVisits = 1)
stavanger50Sites16Vis2 <- weightedDetection(occWeights = occWeights,
visWeights = stavangerVisWeights,
noOccur = 500,
noLocations = 50,
detectProb = 0.426*1.5/16*2,
nIter = 999,
noVisits = 1)
save(stavanger47Sites, stavanger47Sites16, stavanger100Sites16, stavanger66Sites16, stavanger200Sites16, stavanger140Sites16, stavanger50Sites16Vis2, file = "survey_calc_stavanger.Rdata")
surveyFiles <- list.files(pattern = "survey_calc*")
lapply(surveyFiles, load, .GlobalEnv)
surveyFilesToSave <- grep("*Sites*", ls(), value = T)
surveyFilesToSave
save(list = surveyFilesToSave, file = "surveyCalc.Rdata")
data(surveyCalc)
surTab <- tibble("Areal" = rep(c("Norge",
"Oslo",
"Bergen",
"Trondheim",
"Stavanger"), times = 6),
"Antall lokaliteter" = rep(c(47,
66,
100,
140,
200,
50), each = 5),
"Antall besøk" = rep(c(1,
1,
1,
1,
1,
2), each = 5),
"Deteksjonsrate 250x250m" = rep(c(0.426,
0.426 * 1.5,
0.426,
0.426*1.5,
0.426,
0.426*1.5), each = 5))
surTab <- surTab %>%
mutate("Deteksjonerate 1x1km" = `Deteksjonsrate 250x250m` /16 * `Antall besøk`,
"Oppdagbarhet" = c(norge47Sites16$probObs, oslo47Sites16$probObs, bergen47Sites16$probObs, trondheim47Sites16$probObs, stavanger47Sites16$probObs,
norge66Sites16$probObs, oslo66Sites16$probObs, bergen66Sites16$probObs, trondheim66Sites16$probObs, stavanger66Sites16$probObs,
norge100Sites16$probObs, oslo100Sites16$probObs, bergen100Sites16$probObs, trondheim100Sites16$probObs, stavanger100Sites16$probObs,
norge140Sites16$probObs, oslo140Sites16$probObs, bergen140Sites16$probObs, trondheim140Sites16$probObs, stavanger140Sites16$probObs,
norge200Sites16$probObs, oslo200Sites16$probObs, bergen200Sites16$probObs, trondheim200Sites16$probObs, stavanger200Sites16$probObs,
norge50Sites16Vis2$probObs, oslo50Sites16Vis2$probObs, bergen50Sites16Vis2$probObs, trondheim50Sites16Vis2$probObs, stavanger50Sites16Vis2$probObs
))
#surTab
\blandscape
xtable(surTab, caption = "Detection probabilities of a handful of survey regimes, with detection probabilities from an empiric survey.")
\elandscape
How much less is the prob using a subarea? (not much difference it seems)
Calculate the same for 6 year period (cannot split 5 years in two, for the alternative of repeated surveys each year.)
Run calculations for 6 years.
oslo47Sites6Years <- weightedDetection(occWeights = occWeights,
visWeights = osloVisWeights,
noOccur = 500,
noLocations = 47*6,
detectProb = 0.426,
nIter = 999,
noVisits = 1)
oslo47Sites166Years <- weightedDetection(occWeights = occWeights,
visWeights = osloVisWeights,
noOccur = 500,
noLocations = 47*6,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
oslo100Sites166Years <- weightedDetection(occWeights = occWeights,
visWeights = osloVisWeights,
noOccur = 500,
noLocations = 100*6,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
oslo66Sites166Years <- weightedDetection(occWeights = occWeights,
visWeights = osloVisWeights,
noOccur = 500,
noLocations = 66*6,
detectProb = 0.426*1.5/16,
nIter = 999,
noVisits = 1)
oslo200Sites166Years <- weightedDetection(occWeights = occWeights,
visWeights = osloVisWeights,
noOccur = 500,
noLocations = 200*6,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
oslo140Sites166Years <- weightedDetection(occWeights = occWeights,
visWeights = osloVisWeights,
noOccur = 500,
noLocations = 140*6,
detectProb = 0.426*1.5/16,
nIter = 999,
noVisits = 1)
##Bergen
bergen47Sites6Years <- weightedDetection(occWeights = occWeights,
visWeights = bergenVisWeights,
noOccur = 500,
noLocations = 47*6,
detectProb = 0.426,
nIter = 999,
noVisits = 1)
bergen47Sites166Years <- weightedDetection(occWeights = occWeights,
visWeights = bergenVisWeights,
noOccur = 500,
noLocations = 47*6,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
bergen100Sites166Years <- weightedDetection(occWeights = occWeights,
visWeights = bergenVisWeights,
noOccur = 500,
noLocations = 100*6,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
bergen66Sites166Years <- weightedDetection(occWeights = occWeights,
visWeights = bergenVisWeights,
noOccur = 500,
noLocations = 66*6,
detectProb = 0.426*1.5/16,
nIter = 999,
noVisits = 1)
bergen200Sites166Years <- weightedDetection(occWeights = occWeights,
visWeights = bergenVisWeights,
noOccur = 500,
noLocations = 200*6,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
bergen140Sites166Years <- weightedDetection(occWeights = occWeights,
visWeights = bergenVisWeights,
noOccur = 500,
noLocations = 140*6,
detectProb = 0.426*1.5/16,
nIter = 999,
noVisits = 1)
##Trondheim
trondheim47Sites6Years <- weightedDetection(occWeights = occWeights,
visWeights = trondheimVisWeights,
noOccur = 500,
noLocations = 47*6,
detectProb = 0.426,
nIter = 999,
noVisits = 1)
trondheim47Sites166Years <- weightedDetection(occWeights = occWeights,
visWeights = trondheimVisWeights,
noOccur = 500,
noLocations = 47*6,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
trondheim100Sites166Years <- weightedDetection(occWeights = occWeights,
visWeights = trondheimVisWeights,
noOccur = 500,
noLocations = 100*6,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
trondheim66Sites166Years <- weightedDetection(occWeights = occWeights,
visWeights = trondheimVisWeights,
noOccur = 500,
noLocations = 66*6,
detectProb = 0.426*1.5/16,
nIter = 999,
noVisits = 1)
trondheim200Sites166Years <- weightedDetection(occWeights = occWeights,
visWeights = trondheimVisWeights,
noOccur = 500,
noLocations = 200*6,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
trondheim140Sites166Years <- weightedDetection(occWeights = occWeights,
visWeights = trondheimVisWeights,
noOccur = 500,
noLocations = 140*6,
detectProb = 0.426*1.5/16,
nIter = 999,
noVisits = 1)
norge140Sites166Years <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = 140*6,
detectProb = 0.426*1.5/16,
nIter = 999,
noVisits = 1)
#Norge
norge47Sites6Years <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = 47*6,
detectProb = 0.426,
nIter = 999,
noVisits = 1)
norge47Sites166Years <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = 47*6,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
norge100Sites166Years <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = 100*6,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
norge66Sites166Years <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = 66*6,
detectProb = 0.426*1.5/16,
nIter = 999,
noVisits = 1)
norge200Sites166Years <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = 200*6,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
##small expensive
#2 visits, 1.5 detection prob
norge50Sites16Vis26Years <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = 50*6,
detectProb = 0.426*1.5/16*2,
nIter = 999,
noVisits = 1)
oslo50Sites16Vis26Years <- weightedDetection(occWeights = occWeights,
visWeights = osloVisWeights,
noOccur = 500,
noLocations = 50*6,
detectProb = 0.426*1.5/16*2,
nIter = 999,
noVisits = 1)
bergen50Sites16Vis26Years <- weightedDetection(occWeights = occWeights,
visWeights = bergenVisWeights,
noOccur = 500,
noLocations = 50*6,
detectProb = 0.426*1.5/16*2,
nIter = 999,
noVisits = 1)
trondheim50Sites16Vis26Years <- weightedDetection(occWeights = occWeights,
visWeights = trondheimVisWeights,
noOccur = 500,
noLocations = 50*6,
detectProb = 0.426*1.5/16*2,
nIter = 999,
noVisits = 1)
save(oslo47Sites6Years, oslo47Sites166Years, oslo66Sites166Years, oslo100Sites166Years, oslo140Sites166Years, oslo200Sites166Years, oslo50Sites16Vis26Years,
bergen47Sites6Years, bergen47Sites166Years, bergen66Sites166Years, bergen100Sites166Years, bergen140Sites166Years, bergen200Sites166Years, bergen50Sites16Vis26Years,
trondheim47Sites6Years, trondheim47Sites166Years, trondheim66Sites166Years, trondheim100Sites166Years, trondheim140Sites166Years, trondheim200Sites166Years, trondheim50Sites16Vis26Years,
norge47Sites6Years, norge47Sites166Years, norge66Sites166Years, norge100Sites166Years, norge140Sites166Years, norge200Sites166Years, norge50Sites16Vis26Years,
file = "survey_calc6years.Rdata")
stavanger47Sites6Years <- weightedDetection(occWeights = occWeights,
visWeights = stavangerVisWeights,
noOccur = 500,
noLocations = 47*6,
detectProb = 0.426,
nIter = 999,
noVisits = 1)
stavanger47Sites166Years <- weightedDetection(occWeights = occWeights,
visWeights = stavangerVisWeights,
noOccur = 500,
noLocations = 47*6,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
stavanger100Sites166Years <- weightedDetection(occWeights = occWeights,
visWeights = stavangerVisWeights,
noOccur = 500,
noLocations = 100*6,
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
stavanger66Sites166Years <- weightedDetection(occWeights = occWeights,
visWeights = stavangerVisWeights,
noOccur = 500,
noLocations = 66*6,
detectProb = 0.426*1.5/16,
nIter = 999,
noVisits = 1)
#Not enough places with positive occurrence probabilities
# stavanger200Sites166Years <- weightedDetection(occWeights = occWeights,
# visWeights = stavangerVisWeights,
# noOccur = 500,
# noLocations = 200*6,
# detectProb = 0.426/16,
# nIter = 999,
# noVisits = 1)
stavanger200Sites166Years <- tibble("noLocations" = 600, "probObs" = NA)
stavanger140Sites166Years <- weightedDetection(occWeights = occWeights,
visWeights = stavangerVisWeights,
noOccur = 500,
noLocations = 140*6,
detectProb = 0.426*1.5/16,
nIter = 999,
noVisits = 1)
stavanger50Sites16Vis26Years <- weightedDetection(occWeights = occWeights,
visWeights = stavangerVisWeights,
noOccur = 500,
noLocations = 50*6,
detectProb = 0.426*1.5/16*2,
nIter = 999,
noVisits = 1)
save(stavanger47Sites6Years, stavanger47Sites166Years, stavanger66Sites166Years, stavanger100Sites166Years, stavanger140Sites166Years, stavanger200Sites166Years, stavanger50Sites16Vis26Years, file = "survey_calc6yearsStavanger.Rdata")
surTab6Years <- tibble("Areal" = rep(c("Norge",
"Oslo",
"Bergen",
"Trondheim",
"Stavanger"), times = 6),
"Antall lokaliteter per år" = rep(c(47,
66,
100,
140,
200,
50), each = 5),
"Antall lokaliteter totalt" = rep(c(47*6,
66*6,
100*6,
140*6,
200*6,
50*6), each = 5),
"Antall besøk" = rep(c(1,
1,
1,
1,
1,
2), each = 5),
"Deteksjonsrate 250x250m" = rep(c(0.426,
0.426 * 1.5,
0.426,
0.426*1.5,
0.426,
0.426*1.5), each = 5))
surTab6Years <- surTab6Years %>%
mutate("Deteksjonerate 1x1km" = `Deteksjonsrate 250x250m` /16 * `Antall besøk`,
"Oppdagbarhet" = c(norge47Sites166Years$probObs, oslo47Sites166Years$probObs, bergen47Sites166Years$probObs, trondheim47Sites166Years$probObs, stavanger47Sites166Years$probObs,
norge66Sites166Years$probObs, oslo66Sites166Years$probObs, bergen66Sites166Years$probObs, trondheim66Sites166Years$probObs, stavanger66Sites166Years$probObs,
norge100Sites166Years$probObs, oslo100Sites166Years$probObs, bergen100Sites166Years$probObs, trondheim100Sites166Years$probObs, stavanger100Sites166Years$probObs,
norge140Sites166Years$probObs, oslo140Sites166Years$probObs, bergen140Sites166Years$probObs, trondheim140Sites166Years$probObs, stavanger140Sites166Years$probObs,
norge200Sites166Years$probObs, oslo200Sites166Years$probObs, bergen200Sites166Years$probObs, trondheim200Sites166Years$probObs, stavanger200Sites166Years$probObs,
norge50Sites16Vis26Years$probObs, oslo50Sites16Vis26Years$probObs, bergen50Sites16Vis26Years$probObs, trondheim50Sites16Vis26Years$probObs, stavanger50Sites16Vis26Years$probObs
))
#surTab
\blandscape
xtable(surTab6Years, caption = "Detection probabilities of a handful of survey regimes, with detection probabilities from an empiric survey. 6 years of surveying visiting each location only one year.")
\elandscape
Calculate the same but where we revisit each location 2 times. With 47 lokations per year, this gives us in total 47*3 = 141 locations but double the area covered within the same location. We can also measure the spread of known distibutions between the 3 years.
oslo47Sites3Years <- weightedDetection(occWeights = occWeights,
visWeights = osloVisWeights,
noOccur = 500,
noLocations = 47*3,
detectProb = 2*0.426,
nIter = 999,
noVisits = 1)
oslo47Sites163Years <- weightedDetection(occWeights = occWeights,
visWeights = osloVisWeights,
noOccur = 500,
noLocations = 47*3,
detectProb = 2*0.426/16,
nIter = 999,
noVisits = 1)
oslo100Sites163Years <- weightedDetection(occWeights = occWeights,
visWeights = osloVisWeights,
noOccur = 500,
noLocations = 100*3,
detectProb = 2*0.426/16,
nIter = 999,
noVisits = 1)
oslo66Sites163Years <- weightedDetection(occWeights = occWeights,
visWeights = osloVisWeights,
noOccur = 500,
noLocations = 66*3,
detectProb = 2*0.426*1.5/16,
nIter = 999,
noVisits = 1)
oslo200Sites163Years <- weightedDetection(occWeights = occWeights,
visWeights = osloVisWeights,
noOccur = 500,
noLocations = 200*3,
detectProb = 2*0.426/16,
nIter = 999,
noVisits = 1)
oslo140Sites163Years <- weightedDetection(occWeights = occWeights,
visWeights = osloVisWeights,
noOccur = 500,
noLocations = 140*3,
detectProb = 2*0.426*1.5/16,
nIter = 999,
noVisits = 1)
##Bergen
bergen47Sites3Years <- weightedDetection(occWeights = occWeights,
visWeights = bergenVisWeights,
noOccur = 500,
noLocations = 47*3,
detectProb = 2*0.426,
nIter = 999,
noVisits = 1)
bergen47Sites163Years <- weightedDetection(occWeights = occWeights,
visWeights = bergenVisWeights,
noOccur = 500,
noLocations = 47*3,
detectProb = 2*0.426/16,
nIter = 999,
noVisits = 1)
bergen100Sites163Years <- weightedDetection(occWeights = occWeights,
visWeights = bergenVisWeights,
noOccur = 500,
noLocations = 100*3,
detectProb = 2*0.426/16,
nIter = 999,
noVisits = 1)
bergen66Sites163Years <- weightedDetection(occWeights = occWeights,
visWeights = bergenVisWeights,
noOccur = 500,
noLocations = 66*3,
detectProb = 2*0.426*1.5/16,
nIter = 999,
noVisits = 1)
bergen200Sites163Years <- weightedDetection(occWeights = occWeights,
visWeights = bergenVisWeights,
noOccur = 500,
noLocations = 200*3,
detectProb = 2*0.426/16,
nIter = 999,
noVisits = 1)
bergen140Sites163Years <- weightedDetection(occWeights = occWeights,
visWeights = bergenVisWeights,
noOccur = 500,
noLocations = 140*3,
detectProb = 2*0.426*1.5/16,
nIter = 999,
noVisits = 1)
##Trondheim
trondheim47Sites3Years <- weightedDetection(occWeights = occWeights,
visWeights = trondheimVisWeights,
noOccur = 500,
noLocations = 47*3,
detectProb = 2*0.426,
nIter = 999,
noVisits = 1)
trondheim47Sites163Years <- weightedDetection(occWeights = occWeights,
visWeights = trondheimVisWeights,
noOccur = 500,
noLocations = 47*3,
detectProb = 2*0.426/16,
nIter = 999,
noVisits = 1)
trondheim100Sites163Years <- weightedDetection(occWeights = occWeights,
visWeights = trondheimVisWeights,
noOccur = 500,
noLocations = 100*3,
detectProb = 2*0.426/16,
nIter = 999,
noVisits = 1)
trondheim66Sites163Years <- weightedDetection(occWeights = occWeights,
visWeights = trondheimVisWeights,
noOccur = 500,
noLocations = 66*3,
detectProb = 2*0.426*1.5/16,
nIter = 999,
noVisits = 1)
trondheim200Sites163Years <- weightedDetection(occWeights = occWeights,
visWeights = trondheimVisWeights,
noOccur = 500,
noLocations = 200*3,
detectProb = 2*0.426/16,
nIter = 999,
noVisits = 1)
trondheim140Sites163Years <- weightedDetection(occWeights = occWeights,
visWeights = trondheimVisWeights,
noOccur = 500,
noLocations = 140*3,
detectProb = 2*0.426*1.5/16,
nIter = 999,
noVisits = 1)
norge140Sites163Years <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = 140*3,
detectProb = 2*0.426*1.5/16,
nIter = 999,
noVisits = 1)
#Norge
norge47Sites3Years <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = 47*3,
detectProb = 2*0.426,
nIter = 999,
noVisits = 1)
norge47Sites163Years <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = 47*3,
detectProb = 2*0.426/16,
nIter = 999,
noVisits = 1)
norge100Sites163Years <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = 100*3,
detectProb = 2*0.426/16,
nIter = 999,
noVisits = 1)
norge66Sites163Years <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = 66*3,
detectProb = 2*0.426*1.5/16,
nIter = 999,
noVisits = 1)
norge200Sites163Years <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = 200*3,
detectProb = 2*0.426/16,
nIter = 999,
noVisits = 1)
##small expensive
#2 visits, 1.5 detection prob
norge50Sites16Vis23Years <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = 50*3,
detectProb = 2*2*0.426*1.5/16,
nIter = 999,
noVisits = 1)
save(oslo47Sites3Years, oslo47Sites163Years, oslo66Sites163Years, oslo100Sites163Years, oslo140Sites163Years, oslo200Sites163Years, oslo50Sites16Vis23Years,
bergen47Sites3Years, bergen47Sites163Years, bergen66Sites163Years, bergen100Sites163Years, bergen140Sites163Years, bergen200Sites163Years, bergen50Sites16Vis23Years,
trondheim47Sites3Years, trondheim47Sites163Years, trondheim66Sites163Years, trondheim100Sites163Years, trondheim140Sites163Years, trondheim200Sites163Years, trondheim50Sites16Vis23Years,
norge47Sites3Years, norge47Sites163Years, norge66Sites163Years, norge100Sites163Years, norge140Sites163Years, norge200Sites163Years, norge50Sites16Vis23Years,
file = "survey_calc3Years.Rdata")
stavanger47Sites3Years <- weightedDetection(occWeights = occWeights,
visWeights = stavangerVisWeights,
noOccur = 500,
noLocations = 47*3,
detectProb = 2*0.426,
nIter = 999,
noVisits = 1)
stavanger47Sites163Years <- weightedDetection(occWeights = occWeights,
visWeights = stavangerVisWeights,
noOccur = 500,
noLocations = 47*3,
detectProb = 2*0.426/16,
nIter = 999,
noVisits = 1)
stavanger100Sites163Years <- weightedDetection(occWeights =occWeights,
visWeights = stavangerVisWeights,
noOccur = 500,
noLocations = 100*3,
detectProb = 2*0.426/16,
nIter = 999,
noVisits = 1)
stavanger66Sites163Years <- weightedDetection(occWeights = occWeights,
visWeights = stavangerVisWeights,
noOccur = 500,
noLocations = 66*3,
detectProb = 2*0.426*1.5/16,
nIter = 999,
noVisits = 1)
stavanger200Sites163Years <- weightedDetection(occWeights = occWeights,
visWeights = stavangerVisWeights,
noOccur = 500,
noLocations = 200*3,
detectProb = 2*0.426/16,
nIter = 999,
noVisits = 1)
stavanger140Sites163Years <- weightedDetection(occWeights = occWeights,
visWeights = stavangerVisWeights,
noOccur = 500,
noLocations = 140*3,
detectProb = 2*0.426*1.5/16,
nIter = 999,
noVisits = 1)
stavanger50Sites16Vis23Years <- weightedDetection(occWeights = occWeights,
visWeights = stavangerVisWeights,
noOccur = 500,
noLocations = 50*3,
detectProb = 2*2*0.426*1.5/16,
nIter = 999,
noVisits = 1)
save(stavanger47Sites3Years, stavanger47Sites163Years, stavanger100Sites163Years, stavanger66Sites163Years, stavanger200Sites163Years, stavanger140Sites163Years, stavanger50Sites16Vis23Years, file = "survey_calc3YearsStavanger.Rdata")
surTab3years <- tibble("Areal" = rep(c("Norge",
"Oslo",
"Bergen",
"Trondheim",
"Stavanger"), times = 6),
"Antall lokaliteter per år" = rep(c(47,
66,
100,
140,
200,
50), each = 5),
"Antall lokaliteter totalt" = rep(c(47*3,
66*3,
100*3,
140*3,
200*3,
50*3), each = 5),
"Antall besøk" = rep(c(1*2,
1*2,
1*2,
1*2,
1*2,
2*2), each = 5),
"Deteksjonsrate 250x250m" = rep(c(0.426,
0.426 * 1.5,
0.426,
0.426*1.5,
0.426,
0.426*1.5), each = 5))
surTab3years <- surTab3years %>%
mutate("Deteksjonerate 1x1km" = `Deteksjonsrate 250x250m` /16 * `Antall besøk`,
"Oppdagbarhet" = c(norge47Sites163Years$probObs, oslo47Sites163Years$probObs, bergen47Sites163Years$probObs, trondheim47Sites163Years$probObs, stavanger47Sites163Years$probObs,
norge66Sites163Years$probObs, oslo66Sites163Years$probObs, bergen66Sites163Years$probObs, trondheim66Sites163Years$probObs, stavanger66Sites163Years$probObs,
norge100Sites163Years$probObs, oslo100Sites163Years$probObs, bergen100Sites163Years$probObs, trondheim100Sites163Years$probObs, stavanger100Sites163Years$probObs,
norge140Sites163Years$probObs, oslo140Sites163Years$probObs, bergen140Sites163Years$probObs, trondheim140Sites163Years$probObs, stavanger140Sites163Years$probObs,
norge200Sites163Years$probObs, oslo200Sites163Years$probObs, bergen200Sites163Years$probObs, trondheim200Sites163Years$probObs, stavanger200Sites163Years$probObs,
norge50Sites16Vis23Years$probObs, oslo50Sites16Vis23Years$probObs, bergen50Sites16Vis23Years$probObs, trondheim50Sites16Vis23Years$probObs, stavanger50Sites16Vis23Years$probObs
))
#surTab
\blandscape
xtable(surTab3years, caption = "Detection probabilities of a handful of survey regimes, with detection probabilities from an empiric survey. 6 year period visiting each location 2 years.")
\elandscape
Calculate the number of locations and visits needed given the empirical detection rates.
norgeDet0.426Vis1 <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = seq(100, 1000, by = 100),
detectProb = 0.426/16,
nIter = 999,
noVisits = 1)
norgeDet0.426Vis4 <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = seq(100, 1000, by = 100),
detectProb = 0.426/16*4,
nIter = 999,
noVisits = 1)
norgeDet0.426Vis10 <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = seq(100, 1000, by = 100),
detectProb = 0.426/16*10,
nIter = 999,
noVisits = 1)
norgeDet0.639Vis1 <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = seq(100, 1000, by = 100),
detectProb = 0.639/16,
nIter = 999,
noVisits = 1)
norgeDet0.639Vis4 <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = seq(100, 1000, by = 100),
detectProb = 0.639/16*4,
nIter = 999,
noVisits = 1)
norgeDet0.639Vis10 <- weightedDetection(occWeights = occWeights,
visWeights = occWeights,
noOccur = 500,
noLocations = seq(100, 1000, by = 100),
detectProb = 0.639/16*10,
nIter = 999,
noVisits = 1)
save(norgeDet0.426Vis1, norgeDet0.426Vis4, norgeDet0.426Vis10, norgeDet0.639Vis1, norgeDet0.639Vis4, norgeDet0.639Vis10, file = "neededLocEmp.Rdata")
data(neededLocEmp)
plot(norgeDet0.426Vis1, threshold = 0.8,
xlab = "Antall lokaliteter",
ylab = "Observasjonssannsynlighet")
plot(norgeDet0.426Vis4, threshold = 0.8,
xlab = "Antall lokaliteter",
ylab = "Observasjonssannsynlighet")
plot(norgeDet0.426Vis10, threshold = 0.8,
xlab = "Antall lokaliteter",
ylab = "Observasjonssannsynlighet")
plot(norgeDet0.639Vis1, threshold = 0.8,
xlab = "Antall lokaliteter",
ylab = "Observasjonssannsynlighet")
plot(norgeDet0.639Vis4, threshold = 0.8,
xlab = "Antall lokaliteter",
ylab = "Observasjonssannsynlighet")
plot(norgeDet0.639Vis10, threshold = 0.8,
xlab = "Antall lokaliteter",
ylab = "Observasjonssannsynlighet")
\clearpage
References
NINAnor/SurveyPower documentation built on May 7, 2019, 10:43 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
options(list("xtable.include.rownames" = F, "xtable.comment" = F)) require(knitr) knitr::opts_chunk$set(echo = TRUE, tidy = 'styler', dev = c("png", "pdf", "svg"), dev.args=list(family='Arial'), dpi = 600, fig.path = "figure/") require(extrafont) #loadfonts()
Intro
We here explore some effects of rare occurrences and imperfect detection in observation of species, and explore alternative survey regimes. This is mostly relevant for rare species, and not for the general survey of insects. A use case is early detection of alien species, which we explore a little deeper.
The objective here is to identify suitable sampling strategies for a couple of different scenarios, where we maximise the chances of detecting a species at least once while considering the economic costs. The basic principle here is that a species can be present in a sampling location with a given probability of occurrence. The observers then have a given chance (probability) of detecting a species on a visit, if it is present. The resulting probability of detecting a species will thus be dependent on both the probability of occurrence and the probability of detection, how many locations are visited, and how many visits we make at each location.
If $\psi$ represents the occurrence probability and $\theta$ is the detection probability, the probability of observing the species at least once in $J$ locations visited $K$ times is: $$probObserv = 1 - (1 - \psi * (1 - (1 - \theta)^K))^J.$$ To simplify, we here assume that the costs scale roughly linearly to the total amount of visits, i.e. $totalCost = J * K * c$, where $c$ is,the cost of one visit/survey.
First look
We can explore how the overall probability of observing the species, and the associated costs depend on $\psi, \theta, J, K, c$ graphically. I have made a convenience function obsProb to calculate the values. A lot of the results seem pretty obvious in retrospect, but it is still worth plotting them to get a better feel for the possibilities. At the moment, we only consider one and the same cost for each survey visit. It may be worth while coding up using a higher initial cost of the first visit, and lower costs for subsequent visits, to account for the economics of scale.
require(SurveyPower) ##custom library, available through devtools::install_github("NINAnor/SurveyPower") require(dplyr) require(gridExtra) require(ggplot2) require(xtable)
We can specify a single or a range of values for the occurrence probability in each location, the detection probabilities, the number of locations visited, the number of visits per location, and the cost of each visit. To illustrate, we here use occurrence probabilities ranging from 0.01 to 1, two different detection probabilities as 0.5, and 0.8. We specify 30 locations, each visited 4 times, at an individual visit cost of 5000.
obsDf <- obsProb(occProb = seq(0.01, 1, by = 0.05), detectProb = c(0.5, 0.8), locations = 30, visits = 1, visitCost = 5000) obsDf
The probability of observing the species at least once is found in the obsProb column, and the total cost in the column totCost. The function returns an object of a specific class with a custom plotting function (see figure \ref{princPlot}). This makes repeated plottings easier. We can specify a grouping variable to split up the lines.
```r"} plot(obsDf, group = "detectProb", xVar = "occProb", yVar = "obsProb", titleVar = c("locations", "visits"), hline = 0.8)
Figure \ref{princPlot} shows how the overall obervation probability is dependent on both occurrence and detection probability. A threshold of a total observation probability of 0.8 is added for comparison. The threshold is reached in this case when the probability of occurrence in each location is 6%. With 30 locations, the overall observation probability rises quite sharply as a result of increased occurrence probability. With only 4 visits per location, the detection probability limits the overall achievable observation probability.
However, these occurrence probabilities are probably unreasonably high for rare species in the real world. If we assume that we search for a truly rare species with an occurrence probability of only 0.001, we find that reaching overall detectabilities above 80% is challenging (Figure \ref{rareLocPlot}).
```r
rareLocDf <- obsProb(occProb = 0.001,
detectProb = seq(0.4, 0.8, by = 0.2),
locations = seq(30, 1000, by = 20),
visits = seq(4, 16, by = 4),
visitCost = 5000)
rareLocDf
```r"} plot(rareLocDf, group = "visits", xVar = "locations", yVar = "obsProb", hline = 0.8) + facet_grid(occProb ~ detectProb)
In figure \ref{rareLocPlot}, we see that observing a very rare species with a high certainty is difficult even with a very large number of visited locations. In these cases, it doesn't really help to visit each location many times, as the overall probability is limited by the number of locations we visit. Figure \ref{rareVisitPlot} shows the results of maximising the number of visits in fixed, but large number of locations, for a very rare species.
```r
rareVisitDf <- obsProb(occProb = 0.001,
detectProb = seq(0.2, 0.8, by = .2),
locations = 250,
visits = seq(1, 11, by = 5),
visitCost = 5000)
rareVisitDf
```r"} plot(rareVisitDf, group = "detectProb", xVar = "visits", yVar = "obsProb", titleVar = c("occProb", "locations"), hline = 0.8)
Although the overall cost increases linearly with the total number of samples (figure \ref{rareCostPlot}), in cases with very rare species, this doesn't mean that the overall detection probability continues to increase indefinitely (figure \ref{rareVisitPlot}) ```r rareCostDf <- obsProb(occProb = 0.001, detectProb = 0.4, locations = seq(50, 250, by = 50), visits = seq(1, 21, by = 5), visitCost = 5000) rareCostDf
```r"} plot(rareCostDf, group = "locations", xVar = "visits", yVar = "totCost", titleVar = "visitCost")
<!-- ```r"} --> <!-- plot(rareVisitDf, group = "detectProb", --> <!-- xVar = "totCost", --> <!-- yVar = "obsProb", --> <!-- titleVar = c("occProb", "locations")) --> <!-- ``` --> Sampling strategies ============== When we deal with a very rare species, it is little use increasing the number of visits to each location (or to spend money maximizing the detection probability), if we only visit a few locations. In these cases, we have to cover a very large number of locations, just to have a chance of visiting a location with a presence (figure \ref{lowOccurrPlot}). Still, for a very rare species such as displayed in figure \ref{lowOccurrPlot}, with an occurrence probability of 0.01%, reaching an overall observation probability of 80% requires more at least 1650 locations, which would be unfeasible for most survey programs. Alternatively, if detectability is low but presence is relatively high, we should focus on increasing the number of revisits per location, instead of trying to cover many locations (figure \ref{lowDetectPlot}). ```r lowOccurrDf <- obsProb(occProb = 0.001, detectProb = 0.8, locations = seq(50, 2000, by = 50), visits = c(1, seq(5, 20, by = 5)), visitCost = 5000) lowOccurrDf
```r"} plot(lowOccurrDf, group = "visits", xVar = "locations", yVar = "obsProb", titleVar = c("occProb", "detectProb"), hline = 0.8)
```r lowDetectDf <- obsProb(occProb = 0.05, detectProb = 0.2, locations = seq(50, 250, by = 50), visits = c(1, seq(5, 20, by = 5)), visitCost = 5000) lowDetectDf
```r"} plot(lowDetectDf, group = "visits", xVar = "locations", yVar = "obsProb", titleVar = c("occProb", "detectProb"), hline = 0.8)
Plausible values for early detection of alien species. ------------- If we deal with for example early detection of alien species, it is reasonable to assume that we only have to consider rather low occurrence probabilities. We can explore our possibilities of observing a species that occurr in between 0.1 to 1% of all studied locations. For detection probabilities, we choose a range of expected, while hopeful, values. We her set the number of locations to 100 and with two visits, as a reasonable setup to complete within a year. ```r guessDf <- obsProb(occProb = c(0.001, 0.005, 0.01, 0.02, 0.04, 0.05, 0.1), detectProb = c(0.1, seq(0.2, 0.8, by = 0.2)), locations = c(50, seq(100, 300, by = 100)), visits = seq(2, 6, by = 2), visitCost = 5000) guessDf
plot(guessDf, group = "detectProb", xVar = "occProb", yVar = "obsProb", hline = 0.8) + facet_grid(visits ~ locations) + theme(axis.text.x = element_text(size = 5))
```r"} plot(guessDf, group = "detectProb", xVar = "occProb", yVar = "obsProb", hline = 0.8) + facet_grid(visits ~ locations) + scale_y_log10() + scale_x_log10() + theme(axis.text.x = element_text(size = 6)) + ylab("Observasjonssannsynlighet") + xlab("Forekomstsannsynlighet") + guides(color=guide_legend(title="Pdeteksjon"))
Unequally distributed probabilities ================== So far, we have only considered equally distributed probabilities for occurrences, meaning that every grid cell has equal probability of containing an occurrence. In the equation $probObserv = (1 - (1 - \psi)^J) * (1 - (1 - \theta)^K))$, $\psi$ designates the probability that a site we visit contains a certain species. This is the parameter which up until now have been the same for all sites. This is not such a good model for many real world situations. Often we have some knowledge, or at least guesses of where a species is more likely to be present. Therefore, the probability that a specific site that you visit contains a certain species will vary between sites. It will depend both on the probability of a site containing a species, and the probability that you visit that site. We can view the probability of a specific site containing a species as a weighted sampling without replacement. For example, if we know there are a 100 sites containing species x, we can calculate the probability that each site contains species x if we know the probability weights. This probability can be written, following @Efraimidis_2006, as $p_i(k) = \frac{w_i} {\sum_{S_j \in V - S} {w_j}}$. For brevity, however, we will simply designate this probability as $Pr[w_i, n]$, where the weights $w_i$ sums to 1, and $n$ is the number of sites with the species. Similarly, we choose which sites to visit as a sampling without replacement, so that the probability of visiting a site $i$ is $Pr[u_i, v]$, where $u_i$ is the visitation weights, which sums to 1, and $v$ is the number of sites you visit. The probability of visiting a site with an alien species then becomes $\psi_i = Pr[w_i, n] * Pr[u_i, v]$, and the probability of visiting any location inhabited by an alien species is $PoccurrVisit = (1 - \prod_{i}^{}(1 - (Pr[w_i, n] * Pr[u_i, v])))$. For simplicity, we can assume that the probability of detection is equal for each site and visit. The probability of detecting an alien species is then $Pdetect = (1 - \prod_{i}^{}(1 - (Pr[w_i, n] * Pr[u_i, v]))) * (1 - (1 - d)^K$, where d is the detection probability, and K is again the number of visits per site. In the best of worlds, $w_i$ and $u_i$ would match up well, so that we visit the most likely sites to contain a species of interest. The overall probability of detection will decrease as the difference between these variables increase, or in simpler terms we visit the wrong places. Also, the probability will go down the more spread out the probability weights are. In the extreme case, with no spread in these probabilities, there is 100 % certainty that a species will be present in location 1, and 100 % certainty that we will visit just that. In this case, the probabilities are 1 for both occurrence and visitation, so that we are certain we visit a site with a presence. In the other end, there might be no information in the weights, so that the occurrences are randomly spread out, and we visit random sites. In that case, $w_i$ are all the same and $u_i$ are all the same, and we end up with the first equation. In reality, we don't know the true weights that the sites will contain a specific species, and we might choose to visit the wrong sites. For this example calculation, we will assume we know the occurrence weights, and visit them accordingly. In other words, that $w_i = u_i$. If we stipulate the desired detection probability (at e.g 0.8), we can calculate how many sites $v$ we need to visit to be able to detect a species that occurr in $n$ sites with weighted probabilities $w_i$, with a specific certainty. So far, the best estimate for the occurrence weights are the occurrence modelling of alien vascular plants from Olsen et al. 2017. Using this as input, and some simple assumptions, we can calculate the needed number of sites. To get test it out, we can use the 10km scale, which isn't so resource intensive. The 1km scale isn't really interactive since it takes to long time to calculate. ```r require(RPostgreSQL) require(DBI) require(rpostgis) require(SurveyPower) require(tidyverse) data(map10km) con <- dbConnect(RPostgreSQL::PostgreSQL(), host = "gisdata-db.nina.no", dbname = "gisdata", user = "postgjest", password = "gjestpost") #pred <- pgGetRast(con, name = c("hotspot_ias", "bigPred1km"), rast = "rast", bands = 1, boundary = NULL) predQ <- "SELECT ssbid, COALESCE(avg(ST_Value(pred.rast, ST_Centroid(ssb.geom))), 0) pred FROM ssb_data_utm33n.ssb_10km ssb, hotspot_ias.\"evenintbigpred1km\" pred WHERE ST_Intersects(ssb.geom, pred.rast) GROUP BY ssbid" pred <- dbGetQuery(con, predQ) pred <- pred %>% mutate(ssbid = as.character(ssbid)) predMap <- map10km %>% left_join(pred, by = c("ssbid" = "ssbid")) predMap$pred[is.na(predMap$pred)] <- 0 plot(predMap["pred"], border = 0.2)
predWeights10km <- predMap %>% select(sites = ssbid, weights = pred) %>% sf::st_set_geometry(NULL) devtools::use_data(predWeights10km)
data(predWeights10km)
But we start with a situation where the occupied patches are distributed randomly, and we visit the sites randomly. In other words, where the occurrence and visitation weights all are equal.
predWeightsZero<- predWeights10km predWeightsZero$weights <- 1 system.time(predProb50Zero <- weightedDetection(occWeights = predWeightsZero, visWeights = predWeights10km, noOccur = 5, noLocations = seq(50, 200, by = 50), noVisits = 1, detectProb = 1)) predProb50Zero
plot(predProb50Zero, threshold = 0.8)
We can see that the probability of visiting a randomly occupied cell in this case starts from about 0.07 and approaches 0.2 as we increase our number of visited locations from 50 to 200. We can quality check this with a simpler function.
test <- function(noOccur = 50, noLocations = 5, nIter = 999){ prop <- function(noLocations. = noLocations, noOccur. = noOccur){ visited <- sample(1:4057, noLocations., replace = F) #number of 10km cells occupied <- sample(1:4057, noOccur, replace = F) any(visited %in% occupied) } sum(replicate(nIter, prop())/nIter) } test()
We get the same results if the occurrence probabilities is distributed according to informative weights, and only the visitations are random (not shown).
But what happens when we have information about the occurence of the species? In effect, we limit the number of potential sites we visit to a smaller value, which have higher probability of housing the species. We use the prediction map to set the occurrences, and visitation probabilities. We continue with the detection probability set to 1, with just 1 visit per site.
system.time(predProb <- weightedDetection(occWeights = predWeights10km, visWeights = predWeights10km, noOccur = 5, noLocations = seq(50, 200, by = 50), noVisits = 1, detectProb = 1, nIter = 999)) predProb
plot(predProb, threshold = 0.8)
We then reach quite high probabilities for observing the species. This of course depend on the quality of the predictions, i.e. our weights. This prediction map is actually quite informative. From the histogram of weights, we can see that a small number of sites have a proportionally high weight. This limits the realised occurrences quite a bit.
hist(predWeights10km$weights)
Calculations based on 1km map
require(RPostgreSQL) require(DBI) require(rpostgis) require(SurveyPower) require(tidyverse) data(map1km) con <- dbConnect(RPostgreSQL::PostgreSQL(), host = "gisdata-db.nina.no", dbname = "gisdata", user = "postgjest", password = "gjestpost") #pred <- pgGetRast(con, name = c("hotspot_ias", "bigPred1km"), rast = "rast", bands = 1, boundary = NULL) predQ <- "SELECT ssbid, ST_Value(pred.rast, ST_Centroid(ssb.geom)) pred FROM ssb_data_utm33n.ssb_1km ssb, hotspot_ias.\"evenintbigpred1km\" pred WHERE ST_Intersects(ssb.geom, pred.rast) " pred <- dbGetQuery(con, predQ) pred <- pred %>% mutate(ssbid = as.character(ssbid)) predMap1km <- map1km %>% left_join(pred, by = c("ssbid" = "ssbid")) predMap1km$pred[is.na(predMap1km$pred)] <- 0 #plot(predMap1km["pred"]) #Slow
predWeights1km <- predMap1km %>% select(sites = ssbid, weights = pred) %>% sf::st_set_geometry(NULL) devtools::use_data(predWeights1km)
Since we now have about 100 times as many potential sites, we would need to increase the occurrences accordingly. But while 50 out of 4057 10x10km squares sounds reasonable for an "early" detection, multiplying this with 100 yields 5000 locations in a 1x1km grid. This sounds like a lot for an early detection. But it puts the 50 occurrences above into perspective. Surveying 10x10km cells with good detection probability is a tall order.
For the 1x1km analysis, the calculations takes to much time to do on the fly, so we pre-calculate them and load the results.
data("predWeights1km")
It is reasonable to assume that we won't reach a higher detection probability than 0.8 for a single visit, and even that is probably high for a truly novel species. But we can explore the range of occupied sites that are reasonable to handle with such a good detection probability.
data(predOccs)
predOcc500Det0.05 <- weightedDetection(occWeights = predWeights1km, visWeights = predWeights1km, noOccur = 500, noLocations = seq(200, 1000, by = 200), noVisits = 1, detectProb = 0.8/16, nIter = 999) predOcc100Det0.8 <- weightedDetection(occWeights = predWeights1km, visWeights = predWeights1km, noOccur = 100, noLocations = seq(50, 600, by = 50), noVisits = 1, detectProb = 0.8, nIter = 999) predOcc500Det0.8 <- weightedDetection(occWeights = predWeights1km, visWeights = predWeights1km, noOccur = 500, noLocations = seq(50, 300, by = 50), noVisits = 1, detectProb = 0.8, nIter = 999) predOcc200Det0.8 <- weightedDetection(occWeights = predWeights1km, visWeights = predWeights1km, noOccur = 200, noLocations = seq(50, 600, by = 50), noVisits = 1, detectProb = 0.8, nIter = 999) predOcc100Det0.05Vis1 <- weightedDetection(occWeights = predWeights1km, visWeights = predWeights1km, noOccur = 100, noLocations = seq(200, 1000, by = 200), noVisits = 1, detectProb = 0.8/16, nIter = 999) predOcc500Det0.5 <- weightedDetection(occWeights = predWeights1km, visWeights = predWeights1km, noOccur = 500, noLocations = seq(200, 1000, by = 200), noVisits = 1, detectProb = 0.8/16*10, nIter = 999) predOcc500Det0.2 <- weightedDetection(occWeights = predWeights1km, visWeights = predWeights1km, noOccur = 500, noLocations = seq(200, 1000, by = 200), noVisits = 1, detectProb = 0.8/16*4, nIter = 999) save(predOcc500Det0.05, predOcc100Det0.8, predOcc500Det0.8, predOcc200Det0.8, predOcc100Det0.05Vis1, predOcc500Det0.5, predOcc500Det0.2, file = "pred_occs.Rdata")
#collect a bunch of small files into one predFiles <- list.files(pattern = "pred*") lapply(predFiles, load, .GlobalEnv) predObjectsToStore <-grep("predOcc*", ls(), value = T) save(list = predObjectsToStore, file = "predOccs.Rdata")
```r "}
plot(predOcc500Det0.8, threshold = 0.8 , xlab = "Antall lokaliteter", ylab = "Observasjonssannsynlighet")
So for the case of 500 occupied cells, we reach the target observation probability after visiting about 200 sites. In similar fashion, we can explore the case with 200 and 100 occupied cells, respectively. ```r xtable(predOcc200Det0.8)
```r"} plot(predOcc200Det0.8, threshold = 0.8, xlab = "Antall lokaliteter", ylab = "Observasjonssannsynlighet")
```r xtable(predOcc100Det0.8)
```r"} plot(predOcc100Det0.8, threshold = 0.8)
In the case of 200 occupied cells, we reach an overall observation probability of 0.8 after about 500 visited sites with an observation probability of 0.8. With only 100 occupied cells, we don't reach the target of 0.8 even after 600 visited locations.
We can see the effect of a lower detection probability.
Surveying 250x250m squares
-------------
In practice, it can be challenging to survey even a 1x1km grid with any respectable observation probability. We might therefore subdivide the squares in smaller units, with the result that our detection probability decreases from simply not covering the place where the species occupies. If a species just occurs in one out of 16 250x250 squares within a 1x1km square and we visit only one such smaller square, our detection probability drops to 1/16 of the former level. How much the detection probability drops of course depend on the aggregation pattern of the species, i.e. how much of the 1x1km cell it is present in. On a tangent, the number of occurrences we consider to be acceptable within an "early detection" framework depends on the scale of the grid cells considered, and the way that the species are aggregated. 500 out of ca 50 000 1x1km cells constitutes occurrences in 1% of the cells. If we accept a 1% occurrence in the about 800 000 250x250 cells, that amounts to 8000 occurrences. As long as these are not extremy aggregated, this could hardly be seen as an early establishment phase.
Case of 500 1x1km cells occupied, but we survey only 250x250m subsquares
---------------
We start with the case of 500 occurrences spread out in the 50 000 1x1km cells, but when we visit just a 16th of these cells once.
```r
plot(predOcc500Det0.05,
threshold = 0.8,
xlab = "Antall lokaliteter",
ylab = "Observasjonssannsynlighet")
It is clear that these conditions does not let us reach the desired detection probability of 0.8 even with a great number of sampled locations.
We can continue to explore the possibilities with a lower total occurrence, for example when a species is present in 100 out of the 50 000 1x1km grid cells.
plot(predOcc100Det0.05Vis1, threshold = 0.8, xlab = "Antall lokaliteter", ylab = "Observasjonssannsynlighet")
This situation leaves us with very slim chances of detecting the species.
We can mitigate these low numbers by surveying the same 1x1km square multiple times. Note that I here assume that the species is stationary in only one place, and we return to different subplots within the 1x1km square, so that we cover an increasing portion of the sample location. The detectability thus increases linearly. For example, if we have a detection probability of 0.8 in a given location, but we only visit 10 1/16th of that location, we get a detection probability of 0.8/1610 = 0.5. If we would instead visit 2 subsquares of 1/16, we would get 0.8/164 = 0.2
plot(predOcc500Det0.5, threshold = 0.8, xlab = "Antall lokaliteter", ylab = "Observasjonssannsynlighet")
plot(predOcc500Det0.2, threshold = 0.8, xlab = "Antall lokaliteter", ylab = "Observasjonssannsynlighet")
Surveys limited to a specific area
One strategy is to limit the surveyed areas to focus the efforts. We may thereby increase our chances to visit a location that is occupied by a rare species, but we also cannot find species that are not present in the area of interest.
This situation can be modelled using the same function weightedDetection, by setting the visitation probability to zero for the locations outside our focus area. Based on the prediction map, I have selected a set of kommunes around Oslo that seems to encompass the areas with highest predictions.
data(surveyCalc)
Oslo region
occWeights <- predWeights1km osloArea <- c("Vestby", "Svelvik", "Lørenskog", "Frogn", "Ås", "Skedsmo", "Oppegård", "Ski", "Hobøl", "Spydeberg", "Sande", "Oslo", "Bærum", "Asker", "Lier", "Øvre Eiker", "Nedre Eiker", "Drammen", "Hof", "Holmestrand", "Nesodden", "Frogn", "Røyken" ) osloAreaKommuneNr <- c("0211", "0711", "0230", "0215", "0214", "0231", "0217", "0213", "0138", "0123", "0713", "0301", "0219", "0220", "0626", "0624", "0625", "0602", "0714", "0702", "0216", "0215", "0627" ) ssbidToVisit <- map1km %>% filter(KOMMUNENUMMER %in% osloAreaKommuneNr) %>% select(ssbid, kommune)
plot(ssbidToVisit["kommune"], key.width = lcm(4), key.pos = 4, border = 0)
osloVisWeights <- occWeights osloVisWeights$weights[!(osloVisWeights$sites %in% ssbidToVisit$ssbid)] <- 0 summary(osloVisWeights$weights) summary(occWeights$weights)
Bergen
One of several possible demarcations of "Bergen area"
bergenKommuneNr <- c("1265", "1259", "1253", "1260", "1263", "1251", "1411", "1252", "1256", "1266", "1264", "1242", "1247", "1201", "1243", "1246", "1245" ) ssbidBergen <- map1km %>% filter(KOMMUNENUMMER %in% bergenKommuneNr) %>% select(ssbid, kommune)
plot(ssbidBergen["kommune"], key.width = lcm(4), key.pos = 4, border = 0)
bergenVisWeights <- occWeights bergenVisWeights$weights[!(bergenVisWeights$sites %in% ssbidBergen$ssbid)] <- 0 summary(bergenVisWeights$weights) summary(occWeights$weights)
Trondheim
trondheimKommuneNr <- c("1638", "1657", "1636", "1664", "1662", "1601", "1663", "1714", "1653") ssbidTrondheim <- map1km %>% filter(KOMMUNENUMMER %in% trondheimKommuneNr) %>% select(ssbid, kommune)
plot(ssbidTrondheim["kommune"], key.width = lcm(4), key.pos = 4, border = 0)
trondheimVisWeights <- occWeights trondheimVisWeights$weights[!(trondheimVisWeights$sites %in% ssbidTrondheim$ssbid)] <- 0 summary(trondheimVisWeights$weights) summary(occWeights$weights)
Stavanger
occWeights <- predWeights1km stavangerArea <- c("Hå", "Time", "Sandnes", "Randaberg", "Stavanger", "Sola", "Klepp" ) stavangerAreaKommuneNr <- c("1119", "1121", "1102", "1127", "1103", "1124", "1120" ) ssbidToVisit <- map1km %>% filter(KOMMUNENUMMER %in% stavangerAreaKommuneNr) %>% select(ssbid, kommune)
plot(ssbidToVisit["kommune"], key.width = lcm(4), key.pos = 4, border = 0)
stavangerVisWeights <- occWeights stavangerVisWeights$weights[!(stavangerVisWeights$sites %in% ssbidToVisit$ssbid)] <- 0 summary(stavangerVisWeights$weights) sum(stavangerVisWeights$weights >0)
Run calculations
oslo47Sites <- weightedDetection(occWeights = occWeights, visWeights = osloVisWeights, noOccur = 500, noLocations = 47, detectProb = 0.426, nIter = 999, noVisits = 1) oslo47Sites16 <- weightedDetection(occWeights = occWeights, visWeights = osloVisWeights, noOccur = 500, noLocations = 47, detectProb = 0.426/16, nIter = 999, noVisits = 1) oslo100Sites16 <- weightedDetection(occWeights = occWeights, visWeights = osloVisWeights, noOccur = 500, noLocations = 100, detectProb = 0.426/16, nIter = 999, noVisits = 1) oslo66Sites16 <- weightedDetection(occWeights = occWeights, visWeights = osloVisWeights, noOccur = 500, noLocations = 66, detectProb = 0.426*1.5/16, nIter = 999, noVisits = 1) oslo200Sites16 <- weightedDetection(occWeights = occWeights, visWeights = osloVisWeights, noOccur = 500, noLocations = 200, detectProb = 0.426/16, nIter = 999, noVisits = 1) oslo140Sites16 <- weightedDetection(occWeights = occWeights, visWeights = osloVisWeights, noOccur = 500, noLocations = 140, detectProb = 0.426*1.5/16, nIter = 999, noVisits = 1) ##Bergen bergen47Sites <- weightedDetection(occWeights = occWeights, visWeights = bergenVisWeights, noOccur = 500, noLocations = 47, detectProb = 0.426, nIter = 999, noVisits = 1) bergen47Sites16 <- weightedDetection(occWeights = occWeights, visWeights = bergenVisWeights, noOccur = 500, noLocations = 47, detectProb = 0.426/16, nIter = 999, noVisits = 1) bergen100Sites16 <- weightedDetection(occWeights = occWeights, visWeights = bergenVisWeights, noOccur = 500, noLocations = 100, detectProb = 0.426/16, nIter = 999, noVisits = 1) bergen66Sites16 <- weightedDetection(occWeights = occWeights, visWeights = bergenVisWeights, noOccur = 500, noLocations = 66, detectProb = 0.426*1.5/16, nIter = 999, noVisits = 1) bergen200Sites16 <- weightedDetection(occWeights = occWeights, visWeights = bergenVisWeights, noOccur = 500, noLocations = 200, detectProb = 0.426/16, nIter = 999, noVisits = 1) bergen140Sites16 <- weightedDetection(occWeights = occWeights, visWeights = bergenVisWeights, noOccur = 500, noLocations = 140, detectProb = 0.426*1.5/16, nIter = 999, noVisits = 1) ##Trondheim trondheim47Sites <- weightedDetection(occWeights = occWeights, visWeights = trondheimVisWeights, noOccur = 500, noLocations = 47, detectProb = 0.426, nIter = 999, noVisits = 1) trondheim47Sites16 <- weightedDetection(occWeights = occWeights, visWeights = trondheimVisWeights, noOccur = 500, noLocations = 47, detectProb = 0.426/16, nIter = 999, noVisits = 1) trondheim100Sites16 <- weightedDetection(occWeights = occWeights, visWeights = trondheimVisWeights, noOccur = 500, noLocations = 100, detectProb = 0.426/16, nIter = 999, noVisits = 1) trondheim66Sites16 <- weightedDetection(occWeights = occWeights, visWeights = trondheimVisWeights, noOccur = 500, noLocations = 66, detectProb = 0.426*1.5/16, nIter = 999, noVisits = 1) trondheim200Sites16 <- weightedDetection(occWeights = occWeights, visWeights = trondheimVisWeights, noOccur = 500, noLocations = 200, detectProb = 0.426/16, nIter = 999, noVisits = 1) trondheim140Sites16 <- weightedDetection(occWeights = occWeights, visWeights = trondheimVisWeights, noOccur = 500, noLocations = 140, detectProb = 0.426*1.5/16, nIter = 999, noVisits = 1) norge140Sites16 <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = 140, detectProb = 0.426*1.5/16, nIter = 999, noVisits = 1) #Norge norge47Sites <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = 47, detectProb = 0.426, nIter = 999, noVisits = 1) norge47Sites16 <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = 47, detectProb = 0.426/16, nIter = 999, noVisits = 1) norge100Sites16 <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = 100, detectProb = 0.426/16, nIter = 999, noVisits = 1) norge66Sites16 <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = 66, detectProb = 0.426*1.5/16, nIter = 999, noVisits = 1) norge200Sites16 <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = 200, detectProb = 0.426/16, nIter = 999, noVisits = 1) ##small expensive #2 visits, 1.5 detection prob norge50Sites16Vis2 <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = 50, detectProb = 0.426*1.5/16*2, nIter = 999, noVisits = 1) oslo50Sites16Vis2 <- weightedDetection(occWeights = occWeights, visWeights = osloVisWeights, noOccur = 500, noLocations = 50, detectProb = 0.426*1.5/16*2, nIter = 999, noVisits = 1) bergen50Sites16Vis2 <- weightedDetection(occWeights = occWeights, visWeights = bergenVisWeights, noOccur = 500, noLocations = 50, detectProb = 0.426*1.5/16*2, nIter = 999, noVisits = 1) trondheim50Sites16Vis2 <- weightedDetection(occWeights = occWeights, visWeights = trondheimVisWeights, noOccur = 500, noLocations = 50, detectProb = 0.426*1.5/16*2, nIter = 999, noVisits = 1) save(oslo47Sites, oslo47Sites16, oslo66Sites16, oslo100Sites16, oslo140Sites16, oslo200Sites16, oslo50Sites16Vis2, bergen47Sites, bergen47Sites16, bergen66Sites16, bergen100Sites16, bergen140Sites16, bergen200Sites16, bergen50Sites16Vis2, trondheim47Sites, trondheim47Sites16, trondheim66Sites16, trondheim100Sites16, trondheim140Sites16, trondheim200Sites16, trondheim50Sites16Vis2, norge47Sites, norge47Sites16, norge66Sites16, norge100Sites16, norge140Sites16, norge200Sites16, norge50Sites16Vis2, file = "survey_calc.Rdata")
stavanger47Sites <- weightedDetection(occWeights = occWeights, visWeights = stavangerVisWeights, noOccur = 500, noLocations = 47, detectProb = 0.426, nIter = 999, noVisits = 1) stavanger47Sites16 <- weightedDetection(occWeights = occWeights, visWeights = stavangerVisWeights, noOccur = 500, noLocations = 47, detectProb = 0.426/16, nIter = 999, noVisits = 1) stavanger100Sites16 <- weightedDetection(occWeights = occWeights, visWeights = stavangerVisWeights, noOccur = 500, noLocations = 100, detectProb = 0.426/16, nIter = 999, noVisits = 1) stavanger66Sites16 <- weightedDetection(occWeights = occWeights, visWeights = stavangerVisWeights, noOccur = 500, noLocations = 66, detectProb = 0.426*1.5/16, nIter = 999, noVisits = 1) stavanger200Sites16 <- weightedDetection(occWeights = occWeights, visWeights = stavangerVisWeights, noOccur = 500, noLocations = 200, detectProb = 0.426/16, nIter = 999, noVisits = 1) stavanger140Sites16 <- weightedDetection(occWeights = occWeights, visWeights = stavangerVisWeights, noOccur = 500, noLocations = 140, detectProb = 0.426*1.5/16, nIter = 999, noVisits = 1) stavanger50Sites16Vis2 <- weightedDetection(occWeights = occWeights, visWeights = stavangerVisWeights, noOccur = 500, noLocations = 50, detectProb = 0.426*1.5/16*2, nIter = 999, noVisits = 1) save(stavanger47Sites, stavanger47Sites16, stavanger100Sites16, stavanger66Sites16, stavanger200Sites16, stavanger140Sites16, stavanger50Sites16Vis2, file = "survey_calc_stavanger.Rdata")
surveyFiles <- list.files(pattern = "survey_calc*") lapply(surveyFiles, load, .GlobalEnv) surveyFilesToSave <- grep("*Sites*", ls(), value = T) surveyFilesToSave save(list = surveyFilesToSave, file = "surveyCalc.Rdata")
data(surveyCalc)
surTab <- tibble("Areal" = rep(c("Norge", "Oslo", "Bergen", "Trondheim", "Stavanger"), times = 6), "Antall lokaliteter" = rep(c(47, 66, 100, 140, 200, 50), each = 5), "Antall besøk" = rep(c(1, 1, 1, 1, 1, 2), each = 5), "Deteksjonsrate 250x250m" = rep(c(0.426, 0.426 * 1.5, 0.426, 0.426*1.5, 0.426, 0.426*1.5), each = 5)) surTab <- surTab %>% mutate("Deteksjonerate 1x1km" = `Deteksjonsrate 250x250m` /16 * `Antall besøk`, "Oppdagbarhet" = c(norge47Sites16$probObs, oslo47Sites16$probObs, bergen47Sites16$probObs, trondheim47Sites16$probObs, stavanger47Sites16$probObs, norge66Sites16$probObs, oslo66Sites16$probObs, bergen66Sites16$probObs, trondheim66Sites16$probObs, stavanger66Sites16$probObs, norge100Sites16$probObs, oslo100Sites16$probObs, bergen100Sites16$probObs, trondheim100Sites16$probObs, stavanger100Sites16$probObs, norge140Sites16$probObs, oslo140Sites16$probObs, bergen140Sites16$probObs, trondheim140Sites16$probObs, stavanger140Sites16$probObs, norge200Sites16$probObs, oslo200Sites16$probObs, bergen200Sites16$probObs, trondheim200Sites16$probObs, stavanger200Sites16$probObs, norge50Sites16Vis2$probObs, oslo50Sites16Vis2$probObs, bergen50Sites16Vis2$probObs, trondheim50Sites16Vis2$probObs, stavanger50Sites16Vis2$probObs )) #surTab
\blandscape
xtable(surTab, caption = "Detection probabilities of a handful of survey regimes, with detection probabilities from an empiric survey.")
\elandscape
How much less is the prob using a subarea? (not much difference it seems)
Calculate the same for 6 year period (cannot split 5 years in two, for the alternative of repeated surveys each year.)
Run calculations for 6 years.
oslo47Sites6Years <- weightedDetection(occWeights = occWeights, visWeights = osloVisWeights, noOccur = 500, noLocations = 47*6, detectProb = 0.426, nIter = 999, noVisits = 1) oslo47Sites166Years <- weightedDetection(occWeights = occWeights, visWeights = osloVisWeights, noOccur = 500, noLocations = 47*6, detectProb = 0.426/16, nIter = 999, noVisits = 1) oslo100Sites166Years <- weightedDetection(occWeights = occWeights, visWeights = osloVisWeights, noOccur = 500, noLocations = 100*6, detectProb = 0.426/16, nIter = 999, noVisits = 1) oslo66Sites166Years <- weightedDetection(occWeights = occWeights, visWeights = osloVisWeights, noOccur = 500, noLocations = 66*6, detectProb = 0.426*1.5/16, nIter = 999, noVisits = 1) oslo200Sites166Years <- weightedDetection(occWeights = occWeights, visWeights = osloVisWeights, noOccur = 500, noLocations = 200*6, detectProb = 0.426/16, nIter = 999, noVisits = 1) oslo140Sites166Years <- weightedDetection(occWeights = occWeights, visWeights = osloVisWeights, noOccur = 500, noLocations = 140*6, detectProb = 0.426*1.5/16, nIter = 999, noVisits = 1) ##Bergen bergen47Sites6Years <- weightedDetection(occWeights = occWeights, visWeights = bergenVisWeights, noOccur = 500, noLocations = 47*6, detectProb = 0.426, nIter = 999, noVisits = 1) bergen47Sites166Years <- weightedDetection(occWeights = occWeights, visWeights = bergenVisWeights, noOccur = 500, noLocations = 47*6, detectProb = 0.426/16, nIter = 999, noVisits = 1) bergen100Sites166Years <- weightedDetection(occWeights = occWeights, visWeights = bergenVisWeights, noOccur = 500, noLocations = 100*6, detectProb = 0.426/16, nIter = 999, noVisits = 1) bergen66Sites166Years <- weightedDetection(occWeights = occWeights, visWeights = bergenVisWeights, noOccur = 500, noLocations = 66*6, detectProb = 0.426*1.5/16, nIter = 999, noVisits = 1) bergen200Sites166Years <- weightedDetection(occWeights = occWeights, visWeights = bergenVisWeights, noOccur = 500, noLocations = 200*6, detectProb = 0.426/16, nIter = 999, noVisits = 1) bergen140Sites166Years <- weightedDetection(occWeights = occWeights, visWeights = bergenVisWeights, noOccur = 500, noLocations = 140*6, detectProb = 0.426*1.5/16, nIter = 999, noVisits = 1) ##Trondheim trondheim47Sites6Years <- weightedDetection(occWeights = occWeights, visWeights = trondheimVisWeights, noOccur = 500, noLocations = 47*6, detectProb = 0.426, nIter = 999, noVisits = 1) trondheim47Sites166Years <- weightedDetection(occWeights = occWeights, visWeights = trondheimVisWeights, noOccur = 500, noLocations = 47*6, detectProb = 0.426/16, nIter = 999, noVisits = 1) trondheim100Sites166Years <- weightedDetection(occWeights = occWeights, visWeights = trondheimVisWeights, noOccur = 500, noLocations = 100*6, detectProb = 0.426/16, nIter = 999, noVisits = 1) trondheim66Sites166Years <- weightedDetection(occWeights = occWeights, visWeights = trondheimVisWeights, noOccur = 500, noLocations = 66*6, detectProb = 0.426*1.5/16, nIter = 999, noVisits = 1) trondheim200Sites166Years <- weightedDetection(occWeights = occWeights, visWeights = trondheimVisWeights, noOccur = 500, noLocations = 200*6, detectProb = 0.426/16, nIter = 999, noVisits = 1) trondheim140Sites166Years <- weightedDetection(occWeights = occWeights, visWeights = trondheimVisWeights, noOccur = 500, noLocations = 140*6, detectProb = 0.426*1.5/16, nIter = 999, noVisits = 1) norge140Sites166Years <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = 140*6, detectProb = 0.426*1.5/16, nIter = 999, noVisits = 1) #Norge norge47Sites6Years <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = 47*6, detectProb = 0.426, nIter = 999, noVisits = 1) norge47Sites166Years <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = 47*6, detectProb = 0.426/16, nIter = 999, noVisits = 1) norge100Sites166Years <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = 100*6, detectProb = 0.426/16, nIter = 999, noVisits = 1) norge66Sites166Years <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = 66*6, detectProb = 0.426*1.5/16, nIter = 999, noVisits = 1) norge200Sites166Years <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = 200*6, detectProb = 0.426/16, nIter = 999, noVisits = 1) ##small expensive #2 visits, 1.5 detection prob norge50Sites16Vis26Years <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = 50*6, detectProb = 0.426*1.5/16*2, nIter = 999, noVisits = 1) oslo50Sites16Vis26Years <- weightedDetection(occWeights = occWeights, visWeights = osloVisWeights, noOccur = 500, noLocations = 50*6, detectProb = 0.426*1.5/16*2, nIter = 999, noVisits = 1) bergen50Sites16Vis26Years <- weightedDetection(occWeights = occWeights, visWeights = bergenVisWeights, noOccur = 500, noLocations = 50*6, detectProb = 0.426*1.5/16*2, nIter = 999, noVisits = 1) trondheim50Sites16Vis26Years <- weightedDetection(occWeights = occWeights, visWeights = trondheimVisWeights, noOccur = 500, noLocations = 50*6, detectProb = 0.426*1.5/16*2, nIter = 999, noVisits = 1) save(oslo47Sites6Years, oslo47Sites166Years, oslo66Sites166Years, oslo100Sites166Years, oslo140Sites166Years, oslo200Sites166Years, oslo50Sites16Vis26Years, bergen47Sites6Years, bergen47Sites166Years, bergen66Sites166Years, bergen100Sites166Years, bergen140Sites166Years, bergen200Sites166Years, bergen50Sites16Vis26Years, trondheim47Sites6Years, trondheim47Sites166Years, trondheim66Sites166Years, trondheim100Sites166Years, trondheim140Sites166Years, trondheim200Sites166Years, trondheim50Sites16Vis26Years, norge47Sites6Years, norge47Sites166Years, norge66Sites166Years, norge100Sites166Years, norge140Sites166Years, norge200Sites166Years, norge50Sites16Vis26Years, file = "survey_calc6years.Rdata")
stavanger47Sites6Years <- weightedDetection(occWeights = occWeights, visWeights = stavangerVisWeights, noOccur = 500, noLocations = 47*6, detectProb = 0.426, nIter = 999, noVisits = 1) stavanger47Sites166Years <- weightedDetection(occWeights = occWeights, visWeights = stavangerVisWeights, noOccur = 500, noLocations = 47*6, detectProb = 0.426/16, nIter = 999, noVisits = 1) stavanger100Sites166Years <- weightedDetection(occWeights = occWeights, visWeights = stavangerVisWeights, noOccur = 500, noLocations = 100*6, detectProb = 0.426/16, nIter = 999, noVisits = 1) stavanger66Sites166Years <- weightedDetection(occWeights = occWeights, visWeights = stavangerVisWeights, noOccur = 500, noLocations = 66*6, detectProb = 0.426*1.5/16, nIter = 999, noVisits = 1) #Not enough places with positive occurrence probabilities # stavanger200Sites166Years <- weightedDetection(occWeights = occWeights, # visWeights = stavangerVisWeights, # noOccur = 500, # noLocations = 200*6, # detectProb = 0.426/16, # nIter = 999, # noVisits = 1) stavanger200Sites166Years <- tibble("noLocations" = 600, "probObs" = NA) stavanger140Sites166Years <- weightedDetection(occWeights = occWeights, visWeights = stavangerVisWeights, noOccur = 500, noLocations = 140*6, detectProb = 0.426*1.5/16, nIter = 999, noVisits = 1) stavanger50Sites16Vis26Years <- weightedDetection(occWeights = occWeights, visWeights = stavangerVisWeights, noOccur = 500, noLocations = 50*6, detectProb = 0.426*1.5/16*2, nIter = 999, noVisits = 1) save(stavanger47Sites6Years, stavanger47Sites166Years, stavanger66Sites166Years, stavanger100Sites166Years, stavanger140Sites166Years, stavanger200Sites166Years, stavanger50Sites16Vis26Years, file = "survey_calc6yearsStavanger.Rdata")
surTab6Years <- tibble("Areal" = rep(c("Norge", "Oslo", "Bergen", "Trondheim", "Stavanger"), times = 6), "Antall lokaliteter per år" = rep(c(47, 66, 100, 140, 200, 50), each = 5), "Antall lokaliteter totalt" = rep(c(47*6, 66*6, 100*6, 140*6, 200*6, 50*6), each = 5), "Antall besøk" = rep(c(1, 1, 1, 1, 1, 2), each = 5), "Deteksjonsrate 250x250m" = rep(c(0.426, 0.426 * 1.5, 0.426, 0.426*1.5, 0.426, 0.426*1.5), each = 5)) surTab6Years <- surTab6Years %>% mutate("Deteksjonerate 1x1km" = `Deteksjonsrate 250x250m` /16 * `Antall besøk`, "Oppdagbarhet" = c(norge47Sites166Years$probObs, oslo47Sites166Years$probObs, bergen47Sites166Years$probObs, trondheim47Sites166Years$probObs, stavanger47Sites166Years$probObs, norge66Sites166Years$probObs, oslo66Sites166Years$probObs, bergen66Sites166Years$probObs, trondheim66Sites166Years$probObs, stavanger66Sites166Years$probObs, norge100Sites166Years$probObs, oslo100Sites166Years$probObs, bergen100Sites166Years$probObs, trondheim100Sites166Years$probObs, stavanger100Sites166Years$probObs, norge140Sites166Years$probObs, oslo140Sites166Years$probObs, bergen140Sites166Years$probObs, trondheim140Sites166Years$probObs, stavanger140Sites166Years$probObs, norge200Sites166Years$probObs, oslo200Sites166Years$probObs, bergen200Sites166Years$probObs, trondheim200Sites166Years$probObs, stavanger200Sites166Years$probObs, norge50Sites16Vis26Years$probObs, oslo50Sites16Vis26Years$probObs, bergen50Sites16Vis26Years$probObs, trondheim50Sites16Vis26Years$probObs, stavanger50Sites16Vis26Years$probObs )) #surTab
\blandscape
xtable(surTab6Years, caption = "Detection probabilities of a handful of survey regimes, with detection probabilities from an empiric survey. 6 years of surveying visiting each location only one year.")
\elandscape
Calculate the same but where we revisit each location 2 times. With 47 lokations per year, this gives us in total 47*3 = 141 locations but double the area covered within the same location. We can also measure the spread of known distibutions between the 3 years.
oslo47Sites3Years <- weightedDetection(occWeights = occWeights, visWeights = osloVisWeights, noOccur = 500, noLocations = 47*3, detectProb = 2*0.426, nIter = 999, noVisits = 1) oslo47Sites163Years <- weightedDetection(occWeights = occWeights, visWeights = osloVisWeights, noOccur = 500, noLocations = 47*3, detectProb = 2*0.426/16, nIter = 999, noVisits = 1) oslo100Sites163Years <- weightedDetection(occWeights = occWeights, visWeights = osloVisWeights, noOccur = 500, noLocations = 100*3, detectProb = 2*0.426/16, nIter = 999, noVisits = 1) oslo66Sites163Years <- weightedDetection(occWeights = occWeights, visWeights = osloVisWeights, noOccur = 500, noLocations = 66*3, detectProb = 2*0.426*1.5/16, nIter = 999, noVisits = 1) oslo200Sites163Years <- weightedDetection(occWeights = occWeights, visWeights = osloVisWeights, noOccur = 500, noLocations = 200*3, detectProb = 2*0.426/16, nIter = 999, noVisits = 1) oslo140Sites163Years <- weightedDetection(occWeights = occWeights, visWeights = osloVisWeights, noOccur = 500, noLocations = 140*3, detectProb = 2*0.426*1.5/16, nIter = 999, noVisits = 1) ##Bergen bergen47Sites3Years <- weightedDetection(occWeights = occWeights, visWeights = bergenVisWeights, noOccur = 500, noLocations = 47*3, detectProb = 2*0.426, nIter = 999, noVisits = 1) bergen47Sites163Years <- weightedDetection(occWeights = occWeights, visWeights = bergenVisWeights, noOccur = 500, noLocations = 47*3, detectProb = 2*0.426/16, nIter = 999, noVisits = 1) bergen100Sites163Years <- weightedDetection(occWeights = occWeights, visWeights = bergenVisWeights, noOccur = 500, noLocations = 100*3, detectProb = 2*0.426/16, nIter = 999, noVisits = 1) bergen66Sites163Years <- weightedDetection(occWeights = occWeights, visWeights = bergenVisWeights, noOccur = 500, noLocations = 66*3, detectProb = 2*0.426*1.5/16, nIter = 999, noVisits = 1) bergen200Sites163Years <- weightedDetection(occWeights = occWeights, visWeights = bergenVisWeights, noOccur = 500, noLocations = 200*3, detectProb = 2*0.426/16, nIter = 999, noVisits = 1) bergen140Sites163Years <- weightedDetection(occWeights = occWeights, visWeights = bergenVisWeights, noOccur = 500, noLocations = 140*3, detectProb = 2*0.426*1.5/16, nIter = 999, noVisits = 1) ##Trondheim trondheim47Sites3Years <- weightedDetection(occWeights = occWeights, visWeights = trondheimVisWeights, noOccur = 500, noLocations = 47*3, detectProb = 2*0.426, nIter = 999, noVisits = 1) trondheim47Sites163Years <- weightedDetection(occWeights = occWeights, visWeights = trondheimVisWeights, noOccur = 500, noLocations = 47*3, detectProb = 2*0.426/16, nIter = 999, noVisits = 1) trondheim100Sites163Years <- weightedDetection(occWeights = occWeights, visWeights = trondheimVisWeights, noOccur = 500, noLocations = 100*3, detectProb = 2*0.426/16, nIter = 999, noVisits = 1) trondheim66Sites163Years <- weightedDetection(occWeights = occWeights, visWeights = trondheimVisWeights, noOccur = 500, noLocations = 66*3, detectProb = 2*0.426*1.5/16, nIter = 999, noVisits = 1) trondheim200Sites163Years <- weightedDetection(occWeights = occWeights, visWeights = trondheimVisWeights, noOccur = 500, noLocations = 200*3, detectProb = 2*0.426/16, nIter = 999, noVisits = 1) trondheim140Sites163Years <- weightedDetection(occWeights = occWeights, visWeights = trondheimVisWeights, noOccur = 500, noLocations = 140*3, detectProb = 2*0.426*1.5/16, nIter = 999, noVisits = 1) norge140Sites163Years <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = 140*3, detectProb = 2*0.426*1.5/16, nIter = 999, noVisits = 1) #Norge norge47Sites3Years <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = 47*3, detectProb = 2*0.426, nIter = 999, noVisits = 1) norge47Sites163Years <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = 47*3, detectProb = 2*0.426/16, nIter = 999, noVisits = 1) norge100Sites163Years <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = 100*3, detectProb = 2*0.426/16, nIter = 999, noVisits = 1) norge66Sites163Years <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = 66*3, detectProb = 2*0.426*1.5/16, nIter = 999, noVisits = 1) norge200Sites163Years <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = 200*3, detectProb = 2*0.426/16, nIter = 999, noVisits = 1) ##small expensive #2 visits, 1.5 detection prob norge50Sites16Vis23Years <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = 50*3, detectProb = 2*2*0.426*1.5/16, nIter = 999, noVisits = 1) save(oslo47Sites3Years, oslo47Sites163Years, oslo66Sites163Years, oslo100Sites163Years, oslo140Sites163Years, oslo200Sites163Years, oslo50Sites16Vis23Years, bergen47Sites3Years, bergen47Sites163Years, bergen66Sites163Years, bergen100Sites163Years, bergen140Sites163Years, bergen200Sites163Years, bergen50Sites16Vis23Years, trondheim47Sites3Years, trondheim47Sites163Years, trondheim66Sites163Years, trondheim100Sites163Years, trondheim140Sites163Years, trondheim200Sites163Years, trondheim50Sites16Vis23Years, norge47Sites3Years, norge47Sites163Years, norge66Sites163Years, norge100Sites163Years, norge140Sites163Years, norge200Sites163Years, norge50Sites16Vis23Years, file = "survey_calc3Years.Rdata")
stavanger47Sites3Years <- weightedDetection(occWeights = occWeights, visWeights = stavangerVisWeights, noOccur = 500, noLocations = 47*3, detectProb = 2*0.426, nIter = 999, noVisits = 1) stavanger47Sites163Years <- weightedDetection(occWeights = occWeights, visWeights = stavangerVisWeights, noOccur = 500, noLocations = 47*3, detectProb = 2*0.426/16, nIter = 999, noVisits = 1) stavanger100Sites163Years <- weightedDetection(occWeights =occWeights, visWeights = stavangerVisWeights, noOccur = 500, noLocations = 100*3, detectProb = 2*0.426/16, nIter = 999, noVisits = 1) stavanger66Sites163Years <- weightedDetection(occWeights = occWeights, visWeights = stavangerVisWeights, noOccur = 500, noLocations = 66*3, detectProb = 2*0.426*1.5/16, nIter = 999, noVisits = 1) stavanger200Sites163Years <- weightedDetection(occWeights = occWeights, visWeights = stavangerVisWeights, noOccur = 500, noLocations = 200*3, detectProb = 2*0.426/16, nIter = 999, noVisits = 1) stavanger140Sites163Years <- weightedDetection(occWeights = occWeights, visWeights = stavangerVisWeights, noOccur = 500, noLocations = 140*3, detectProb = 2*0.426*1.5/16, nIter = 999, noVisits = 1) stavanger50Sites16Vis23Years <- weightedDetection(occWeights = occWeights, visWeights = stavangerVisWeights, noOccur = 500, noLocations = 50*3, detectProb = 2*2*0.426*1.5/16, nIter = 999, noVisits = 1) save(stavanger47Sites3Years, stavanger47Sites163Years, stavanger100Sites163Years, stavanger66Sites163Years, stavanger200Sites163Years, stavanger140Sites163Years, stavanger50Sites16Vis23Years, file = "survey_calc3YearsStavanger.Rdata")
surTab3years <- tibble("Areal" = rep(c("Norge", "Oslo", "Bergen", "Trondheim", "Stavanger"), times = 6), "Antall lokaliteter per år" = rep(c(47, 66, 100, 140, 200, 50), each = 5), "Antall lokaliteter totalt" = rep(c(47*3, 66*3, 100*3, 140*3, 200*3, 50*3), each = 5), "Antall besøk" = rep(c(1*2, 1*2, 1*2, 1*2, 1*2, 2*2), each = 5), "Deteksjonsrate 250x250m" = rep(c(0.426, 0.426 * 1.5, 0.426, 0.426*1.5, 0.426, 0.426*1.5), each = 5)) surTab3years <- surTab3years %>% mutate("Deteksjonerate 1x1km" = `Deteksjonsrate 250x250m` /16 * `Antall besøk`, "Oppdagbarhet" = c(norge47Sites163Years$probObs, oslo47Sites163Years$probObs, bergen47Sites163Years$probObs, trondheim47Sites163Years$probObs, stavanger47Sites163Years$probObs, norge66Sites163Years$probObs, oslo66Sites163Years$probObs, bergen66Sites163Years$probObs, trondheim66Sites163Years$probObs, stavanger66Sites163Years$probObs, norge100Sites163Years$probObs, oslo100Sites163Years$probObs, bergen100Sites163Years$probObs, trondheim100Sites163Years$probObs, stavanger100Sites163Years$probObs, norge140Sites163Years$probObs, oslo140Sites163Years$probObs, bergen140Sites163Years$probObs, trondheim140Sites163Years$probObs, stavanger140Sites163Years$probObs, norge200Sites163Years$probObs, oslo200Sites163Years$probObs, bergen200Sites163Years$probObs, trondheim200Sites163Years$probObs, stavanger200Sites163Years$probObs, norge50Sites16Vis23Years$probObs, oslo50Sites16Vis23Years$probObs, bergen50Sites16Vis23Years$probObs, trondheim50Sites16Vis23Years$probObs, stavanger50Sites16Vis23Years$probObs )) #surTab
\blandscape
xtable(surTab3years, caption = "Detection probabilities of a handful of survey regimes, with detection probabilities from an empiric survey. 6 year period visiting each location 2 years.")
\elandscape
Calculate the number of locations and visits needed given the empirical detection rates.
norgeDet0.426Vis1 <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = seq(100, 1000, by = 100), detectProb = 0.426/16, nIter = 999, noVisits = 1) norgeDet0.426Vis4 <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = seq(100, 1000, by = 100), detectProb = 0.426/16*4, nIter = 999, noVisits = 1) norgeDet0.426Vis10 <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = seq(100, 1000, by = 100), detectProb = 0.426/16*10, nIter = 999, noVisits = 1) norgeDet0.639Vis1 <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = seq(100, 1000, by = 100), detectProb = 0.639/16, nIter = 999, noVisits = 1) norgeDet0.639Vis4 <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = seq(100, 1000, by = 100), detectProb = 0.639/16*4, nIter = 999, noVisits = 1) norgeDet0.639Vis10 <- weightedDetection(occWeights = occWeights, visWeights = occWeights, noOccur = 500, noLocations = seq(100, 1000, by = 100), detectProb = 0.639/16*10, nIter = 999, noVisits = 1) save(norgeDet0.426Vis1, norgeDet0.426Vis4, norgeDet0.426Vis10, norgeDet0.639Vis1, norgeDet0.639Vis4, norgeDet0.639Vis10, file = "neededLocEmp.Rdata")
data(neededLocEmp)
plot(norgeDet0.426Vis1, threshold = 0.8, xlab = "Antall lokaliteter", ylab = "Observasjonssannsynlighet")
plot(norgeDet0.426Vis4, threshold = 0.8, xlab = "Antall lokaliteter", ylab = "Observasjonssannsynlighet")
plot(norgeDet0.426Vis10, threshold = 0.8, xlab = "Antall lokaliteter", ylab = "Observasjonssannsynlighet")
plot(norgeDet0.639Vis1, threshold = 0.8, xlab = "Antall lokaliteter", ylab = "Observasjonssannsynlighet")
plot(norgeDet0.639Vis4, threshold = 0.8, xlab = "Antall lokaliteter", ylab = "Observasjonssannsynlighet")
plot(norgeDet0.639Vis10, threshold = 0.8, xlab = "Antall lokaliteter", ylab = "Observasjonssannsynlighet")
\clearpage
References
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.