-
Notifications
You must be signed in to change notification settings - Fork 1
/
Sim_BYrange_YearLast.R
96 lines (80 loc) · 3.53 KB
/
Sim_BYrange_YearLast.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
#' @title Simulate birth year ranges (min-max)
#'
#' @description Simulate birth year minimum -- maximum with varying widths based
#' on an 'actual' birth year.
#'
#' @param LifeHistData dataframe with columns ID - (Sex) - BirthYear. Column
#' name must be 'BirthYear', case sensitive.
#' @param pwidth named vector, where names are widths (whole numbers) and values
#' are probabilities. Passed to \code{sample}.
#' @param delBY logical, delete all exact birth years, i.e. set column
#' \code{BirthYear} to \code{NA}.
#'
#' @return dataframe \code{LifeHistData} with added columns \code{BY.min} and
#' \code{BY.max}.
#'
#' @examples
#' LH_griffin_X <- SimBYR(LH_griffin)
SimBYR <- function(LifeHistData = NULL,
pwidth = c('1'=0.2, '2'=0.5, '3'=0.4, '4'=0.1),
delBY = TRUE)
# pout = 0.01) # probability that true BY is outside range (-1 or +1)
{
# width
W <- as.numeric(
sample(names(pwidth), size = nrow(LifeHistData), prob = pwidth, replace = TRUE) )
# offset
Z <- sapply(seq_along(W), function(i) sample.int(W[i] +1, size=1) -1)
LifeHistData$BY.min <- LifeHistData$BirthYear - Z
LifeHistData$BY.max <- LifeHistData$BY.min + W
if (delBY) LifeHistData$BirthYear <- NA
return( LifeHistData )
}
################################################################################
################################################################################
################################################################################
#' @title Simulate last year of reproduction
#'
#' @description Simulate last year of reproduction, based on actual last year
#' according to pedigree and based on age for non-breeders.
#'
#' @param LifeHistData dataframe with columns ID - Sex - BirthYear. Column
#' names must be exact, case sensitive.
#' @param Pedigree dataframe with columns id - dam - sire. Column
#' names must be exact, case sensitive. Use \code{PedPolish} first if the
#' pedigree is not generated by \code{sequoia}.
#' @param pAfter named vector, where names are years after actual last year
#' (whole numbers), and values are probabilities. 'X' = unknown last year.
#' @param pAge named vector, for individuals with no offspring in
#' \code{Pedigree}: names are age in years, values are probabilities, 'X' =
#' unknown last year.
#'
#' @return dataframe with columns 'ID' and 'YearLast'.
#'
#' @examples
#' LH_griffin_Y <- merge(LH_griffin,
#' SimLY(LH_griffin, Ped_griffin))
#'
SimLY <- function(LifeHistData = NULL,
Pedigree = NULL,
pAfter = c('0'=0.4, '1'=0.3, 'X'=0.3), # have offspring
pAge = c('0'=0.1, '1'=0.3, '2'=0.4, 'X'=0.2)) # no offspring
{
PedY <- merge(Pedigree, LifeHistData, by.x='id', by.y='ID', all.x=TRUE)
LY.true <- NA
for (i in 1:nrow(PedY)) {
BY.off <- PedY$BirthYear[PedY$dam %in% PedY$id[i] | PedY$sire %in% PedY$id[i]]
LY.true[i] <- ifelse(length(BY.off)>0, max(BY.off, na.rm=TRUE), NA)
}
LY.sim <- NA
Xafter <- suppressWarnings( as.integer(names(pAfter)) )
Xage <- suppressWarnings( as.integer(names(pAge)) )
WithOff <- !is.na(LY.true)
NoOff <- is.na(LY.true)
LY.sim[ WithOff ] <- LY.true[ WithOff ] + sample(Xafter, size=sum(WithOff),
prob=pAfter, replace=TRUE)
LY.sim[ NoOff ] <- PedY$BirthYear[ NoOff ] + sample(Xage, size=sum(NoOff),
prob=pAge, replace=TRUE)
return( data.frame(ID = PedY$id,
YearLast = LY.sim) )
}