-
Notifications
You must be signed in to change notification settings - Fork 2
/
rbeta_dmp.R
56 lines (45 loc) · 1.32 KB
/
rbeta_dmp.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
#' Customised sampling functions for the Beta distributions
#'
#' generate random samples from a beta distribution, parameterized as mean and sd, and returning NAs if conditions are not met
#'
#' @param n An integer value. The number of samples to generate
#' @param p A decimal value. The value used to calculate parameters for the beta distribution
#' @param sd A decimal value. The standard deviation of the beta distribution to simulate
#' @return a vector of samples values from the beta distribution
#'
#' @examples
#' rbeta_dmp(n=100,p=0.9,sd=0.01)
#' @export
#'
rbeta_dmp <- function(n, p, sd){
eta <- p*(1-p)/sd^2 - 1
alpha <- eta*p
beta <- eta*(1-p)
betaMeanVarCond <- sd^2 < p*(1-p)
if(is.na(p) | is.na(sd)){
out <- rep(NA, n)
warning("NA values for p and/or sd - NAs produced")
}else{
if(p >= 0 & p <= 1){
if(sd < 0){
out <- rep(NA, n)
warning("stdev < 0 - NAs produced")
}
if(sd == 0){
out <- rep(p, n)
}
if(sd > 0){
if(betaMeanVarCond){
out <- stats::rbeta(n, shape1 = alpha, shape2 = beta)
}else{
out <- rep(NA, n)
warning("condition var < p*(1 - p) not met - NAs produced")
}
}
}else{
out <- rep(NA, n)
warning("p < 0 | p > 1 - NAs produced")
}
}
return(out)
}