Skip to content

Commit

Permalink
Browse files Browse the repository at this point in the history
  • Loading branch information
see24 committed Dec 18, 2023
2 parents 43b1ef2 + 50ca8ee commit c6b380b
Show file tree
Hide file tree
Showing 6 changed files with 28 additions and 18 deletions.
2 changes: 1 addition & 1 deletion R/caribouBayesianIPM.R
Original file line number Diff line number Diff line change
Expand Up @@ -473,7 +473,7 @@ caribouBayesianIPM <- function(survData = system.file("extdata/simSurvData.csv",

sp.params <- c("S.annual.KM", "R", "Rfemale", "pop.growth",
"fpop.size", "var.R.real", "l.R", "l.Saf",
"beta.Rec.anthro", "beta.Rec.fire", "beta.Saf")
"beta.Rec.anthro", "beta.Rec.fire", "beta.Saf","composition.bias")



Expand Down
6 changes: 4 additions & 2 deletions R/getScenarioDefaults.R
Original file line number Diff line number Diff line change
Expand Up @@ -17,6 +17,7 @@
#' @param sQuantile number in 0, 1. Survival quantile
#' @param projYears Number of years of projections
#' @param obsYears Number of years of observations
#' @param preYears Number of years before monitoring begins
#' @param qMin number in 0, 1. Minimum ratio of bulls to cows in composition survey groups.
#' @param qMax number in 0, 1. Maximum ratio of bulls to cows in composition survey groups.
#' @param uMin number in 0, 1. Minimum probability of misidentifying young bulls as adult females and vice versa in composition survey.
Expand Down Expand Up @@ -49,7 +50,8 @@
getScenarioDefaults <- function(paramTable = NULL,
iFire = 0, iAnthro = 0, obsAnthroSlope = 2, projAnthroSlope = 2,
rSlopeMod = 1, sSlopeMod = 1,
rQuantile = 0.5, sQuantile = 0.5, projYears = 35, obsYears = 15, N0 = 1000,
rQuantile = 0.5, sQuantile = 0.5, projYears = 35,
obsYears = 15, preYears=0, N0 = 1000,
adjustR = TRUE, assessmentYrs = 3,qMin=0,qMax =0.6,
uMin = 0, uMax = 0.2, zMin = 0, zMax = 0.2, cowMult = 6,
collarInterval = NA, cowCount = NA,
Expand Down Expand Up @@ -90,7 +92,7 @@ getScenarioDefaults <- function(paramTable = NULL,
}

if (!is.element("startYear", names(paramTable))) {
paramTable$startYear <- paramTable$curYear - paramTable$obsYears + 1
paramTable$startYear <- paramTable$curYear - paramTable$obsYears - paramTable$preYears + 1
}

return(paramTable)
Expand Down
2 changes: 1 addition & 1 deletion R/internal_demog_fns.R
Original file line number Diff line number Diff line change
Expand Up @@ -202,7 +202,7 @@ simCalfCowRatios <- function(cowCounts, minYr, exData) {
ageRatioOut <- merge(ageRatioOut,
subset(exData, select = c("Year", "recruitment","N")))

if(sum(ageRatioOut$cow>ageRatioOut$N)){
if(sum(ageRatioOut$cow>ageRatioOut$N,na.rm=T)){
warning("The expected number of cows in composition survey exceeds population size. Adjusting cows in survey for consistency.")
}

Expand Down
32 changes: 19 additions & 13 deletions R/simulateObservations.R
Original file line number Diff line number Diff line change
Expand Up @@ -90,11 +90,11 @@ simulateObservations <- function(paramTable, cowCounts = NULL,

if(!is.null(cowCounts)){
testTable(cowCounts, c("Year", "Count", "Class"),
req_vals = list(Year = paramTable$startYear:(paramTable$startYear+paramTable$obsYears-1)),
req_vals = list(Year = (paramTable$startYear+paramTable$preYears):(paramTable$startYear+paramTable$preYears+paramTable$obsYears-1)),
acc_vals = list(Class = "cow"))
} else if(!is.null(paramTable$cowCount)){
cowCounts <- data.frame(Year = paramTable$startYear:
(paramTable$startYear+paramTable$obsYears-1),
cowCounts <- data.frame(Year = (paramTable$startYear+paramTable$preYears):
(paramTable$startYear+paramTable$preYears+paramTable$obsYears-1),
Count = paramTable$cowCount,
Class = "cow")
} else if(is.null(paramTable$cowCount) & is.null(paramTable$cowMult)){
Expand All @@ -104,10 +104,10 @@ simulateObservations <- function(paramTable, cowCounts = NULL,

if(!is.null(freqStartsByYear)){
testTable(freqStartsByYear, c("Year", "numStarts"),
acc_vals = list(Year = paramTable$startYear:(paramTable$startYear+paramTable$obsYears-1)))
acc_vals = list(Year = (paramTable$startYear+paramTable$preYears):(paramTable$startYear+paramTable$preYears+paramTable$obsYears-1)))
} else if(!is.null(paramTable$collarCount)){
freqStartsByYear <- data.frame(Year = paramTable$startYear:
(paramTable$startYear+paramTable$obsYears-1),
freqStartsByYear <- data.frame(Year = (paramTable$startYear+paramTable$preYears):
(paramTable$startYear+paramTable$preYears+paramTable$obsYears-1),
numStarts = paramTable$collarCount)
}else {
stop("One of freqStartsByYear or paramTable$collarCount must be provided",
Expand All @@ -117,9 +117,9 @@ simulateObservations <- function(paramTable, cowCounts = NULL,
# Simulate covariate table
if (is.null(distScen)) {
covariates <- simCovariates(paramTable$iAnthro, paramTable$iFire,
paramTable$obsYears + paramTable$projYears,
paramTable$preYears+paramTable$obsYears + paramTable$projYears,
paramTable$obsAnthroSlope, paramTable$projAnthroSlope,
paramTable$obsYears + 1)
paramTable$obsYears + paramTable$preYears + 1)
simDisturbance <- covariates
simDisturbance$Year <- paramTable$startYear + simDisturbance$time - 1

Expand All @@ -132,14 +132,14 @@ simulateObservations <- function(paramTable, cowCounts = NULL,
} else {
simDisturbance <- distScen
simDisturbance$time <- simDisturbance$Year - paramTable$startYear + 1
simDisturbance <- filter(simDisturbance, .data$Year <= (paramTable$startYear + paramTable$obsYears - 1 + paramTable$projYears) &
simDisturbance <- filter(simDisturbance, .data$Year <= (paramTable$startYear + paramTable$preYears + paramTable$obsYears - 1 + paramTable$projYears) &
.data$Year >= paramTable$startYear)
}

# simulate true population trajectory
suppressMessages(
popMetrics <- simTrajectory(
numYears = paramTable$obsYears + paramTable$projYears,
numYears = paramTable$preYears + paramTable$obsYears + paramTable$projYears,
covariates = simDisturbance,
popGrowthTable = populationGrowthTable,
survivalModelNumber = survivalModelNumber,
Expand Down Expand Up @@ -171,15 +171,17 @@ simulateObservations <- function(paramTable, cowCounts = NULL,
names_from = "MetricTypeID",
values_from = "Amount")
popMetricsWide$Year <- paramTable$startYear + popMetricsWide$Timestep - 1

popMetricsWide <- subset(popMetricsWide,popMetricsWide$Timestep > paramTable$preYears)

exData <- subset(popMetricsWide, (popMetricsWide$Timestep <= paramTable$obsYears))
exData <- subset(popMetricsWide, popMetricsWide$Timestep <= (paramTable$obsYears+paramTable$preYears))

# Now apply observation process model to get simulated calf:cow and survival data.
# Use sample sizes in example input data e.g. Eaker

# reduce sim data tables to length of observations prior to max year
minYr <- paramTable$startYear
maxYr <- paramTable$startYear + paramTable$obsYears + paramTable$projYears - 1
minYr <- paramTable$startYear+paramTable$preYears
maxYr <- paramTable$startYear + paramTable$preYears + paramTable$obsYears + paramTable$projYears - 1

# simulate survival data from survival probability.
if (is.element("collarInterval", names(paramTable)) & is.null(freqStartsByYearIn)) {
Expand Down Expand Up @@ -230,6 +232,10 @@ simulateObservations <- function(paramTable, cowCounts = NULL,
file.path(writeFilesDir, paste0("simSurvData", paramTable$label, ".csv")),
row.names = FALSE)
}

#cut simDisturbance table to start at first obs year
simDisturbance = subset(simDisturbance,Year>=minYr)

return(list(minYr = minYr, maxYr = maxYr, simDisturbance = simDisturbance,
simSurvObs = simSurvObs, ageRatioOut = ageRatioOut,
exData = popMetricsWide, paramTable = paramTable))
Expand Down
1 change: 0 additions & 1 deletion inst/templates/JAGS_template.txt
Original file line number Diff line number Diff line change
Expand Up @@ -51,7 +51,6 @@ for(k in 2:nYears){

#NOTE truncate distribution to ensure intercept between 0 and 0.99. See thinkingAboutPriors.R
l.Saf ~ dnorm(l.Saf.Prior1,pow(l.Saf.Prior2,-2)) T(-10,-0.01)
#NOTE: see addInterannualVar fn in caribouMetrics for logic. But note shape of interannual variation distribution is different because using normally distributed random effect.
sig.Saf~dnorm(exp(l.Saf.Prior1)*sig.Saf.Prior1,pow(sig.Saf.Prior2,-2))
tau.Saf<-pow(sig.Saf, -2)
var.Saf<-1/tau.Saf
Expand Down
3 changes: 3 additions & 0 deletions man/getScenarioDefaults.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

0 comments on commit c6b380b

Please sign in to comment.