20

我正在使用 rnorm 模拟数据,但我需要设置上限和下限,有人知道该怎么做吗?

代码:

rnorm(n = 10, mean = 39.74, sd = 25.09)

上限需要340,下限0

我问这个问题是因为我将 SAS 代码重写为 R 代码。我从来没有使用过SAS。我正在尝试重写以下代码:

sim_sample(simtot=100000,seed=10004,lbound=0,ubound=340,round_y=0.01,round_m=0.01,round_sd=0.01,n=15,m=39.74,sd=25.11,mk=4)
4

6 回答 6

30

rtruncnorm() 函数将返回您需要的结果。

  library(truncnorm)
  rtruncnorm(n=10, a=0, b=340, mean=39.4, sd=25.09)
于 2013-10-13T21:34:48.447 回答
13

您可以制作自己的截断法线采样器,不需要您非常简单地丢弃观察结果

rtnorm <- function(n, mean, sd, a = -Inf, b = Inf){
    qnorm(runif(n, pnorm(a, mean, sd), pnorm(b, mean, sd)), mean, sd)
}
于 2013-10-13T21:46:31.377 回答
10

像这样?

mysamp <- function(n, m, s, lwr, upr, nnorm) {
  samp <- rnorm(nnorm, m, s)
  samp <- samp[samp >= lwr & samp <= upr]
  if (length(samp) >= n) {
    return(sample(samp, n))
  }  
  stop(simpleError("Not enough values to sample from. Try increasing nnorm."))
}

set.seed(42)
mysamp(n=10, m=39.74, s=25.09, lwr=0, upr=340, nnorm=1000)
#[1] 58.90437 38.72318 19.64453 20.24153 39.41130 12.80199 59.88558 30.88578 19.66092 32.46025

但是,结果不是正态分布的,通常不会有您指定的平均值和标准差(特别是如果限制围绕指定的平均值不对称)。

编辑:

根据您的评论,您似乎想翻译此 SAS 功能。我不是 SAS 用户,但这应该或多或少相同:

mysamp <- function(n, m, s, lwr, upr, rounding) {
  samp <- round(rnorm(n, m, s), rounding)
  samp[samp < lwr] <- lwr
  samp[samp > upr] <- upr
  samp
}

set.seed(8)
mysamp(n=10, m=39.74, s=25.09, lwr=0, upr=340, rounding=3)
#[1] 37.618 60.826 28.111 25.920 58.207 37.033 35.467 12.434  0.000 24.857

然后,您可能希望使用它replicate来运行模拟。或者,如果您想要更快的代码:

sim <- matrix(mysamp(n=10*10, m=39.74, s=25.09, lwr=0, upr=340, rounding=3), 10)
means <- colMeans(sim)
sds <- apply(sim, 2, sd)
于 2013-10-13T08:52:19.653 回答
0

这是我为实现相同目的而编写的函数。它将函数的结果归一化rnorm,然后对其进行调整以适应范围。

注意:标准偏差和平均值(如果指定)在标准化过程中发生变化。

#' Creates a random normal distribution within the specified bounds.
#' 
#' WARNING: This function does not preserve the standard deviation or mean.
#' @param n The number of values to be generated
#' @param mean The mean of the distribution
#' @param sd The standard deviation of the distribution
#' @param lower The lower limit of the distribution
#' @param upper The upper limit of the distribution
rtnorm <- function(n, mean=NA, sd=1, lower=-1, upper=1){
  mean = ifelse(is.na(mean)|| mean < lower || mean > upper,
                mean(c(lower, upper)), mean)
  data <- rnorm(n, mean=m, sd=sd) # data

  if (!is.na(lower) && !is.na(upper)){ # adjust data to specified range
    drange <- range(data)           # data range
    irange <- range(lower, upper)   # input range
    data <- (data - drange[1])/(drange[2] - drange[1]) # normalize data (make it 0 to 1)
    data <- (data * (irange[2] - irange[1]))+irange[1] # adjust to specified range
  }
  return(data)
}
于 2015-02-22T18:45:31.590 回答
0

假设您想要正好 10 个数字,而不是 >0、<340 的子集(并且夜晚不是正态分布):

    aa <- rnorm(n = 10, mean = 39.74, s = 25.09)

    while(any(aa<0 | aa>340)) { aa <- rnorm(n = 10, mean = 39.74, s = 25.09) }
于 2013-10-13T08:55:19.603 回答
0

有几种方法可以设置正态分布的上限和下限,什么会导致结果不再正态分布

假设 a mean=0sd=1产生N=1e5的值的下边界为LO=-1,上边界为UP=2

N <- 1e5L
LO <- -1
UP <- 2

将异常值移动到边界(@Roland)

set.seed(42)
x <- pmax(LO, pmin(UP, rnorm(N)))
mean(x)
#[1] 0.07238029
median(x)
#[1] -0.002066374
sd(x)
#[1] 0.8457605
hist(x, 30)

将异常值设置为边界的历史记录

剪切 (@Dason, @Roland, truncnorm::rtruncnorm, MCMCglmm::rtnorm) 的异常值

set.seed(42)
x <- qnorm(runif(N, pnorm(LO), pnorm(UP)))
mean(x)
#[1] 0.2317875
median(x)
#[1] 0.173679
sd(x)
#[1] 0.7236536

hist of Cut 异常值

规模(@Alex Essilfie)

set.seed(42)
x <- rnorm(N)
x <- (x-min(x))/(max(x)-min(x))*(UP-LO)+LO
mean(x)
#[1] 0.4474876
median(x)
#[1] 0.4482257
sd(x)
#[1] 0.3595199

规模历史

方法的组合。例如切割和缩放:

set.seed(42)
x <- qnorm(runif(N, pnorm(-3), pnorm(3)))
x <- (x-min(x))/(max(x)-min(x))*(UP-LO)+LO
mean(x)
#[1] 0.5010759
median(x)
#[1] 0.5014713
sd(x)
#[1] 0.4957751

组合历史

不对称组合

set.seed(42)
n <- round(N*abs(LO)/diff(range(c(LO, UP))))
x <- c(qnorm(runif(n, pnorm(-3), 0.5)), qnorm(runif(N-n, 0.5, pnorm(3))))
x <- ifelse(x < 0, x/min(x)*LO, x/max(x)*UP)
mean(x)
#[1] 0.2651627
median(x)
#[1] 0.2127903
sd(x)
#[1] 0.5078264

不对称组合的历史

于 2021-04-07T10:02:52.243 回答