生成具有固定均值和标准差的随机数

Question

当使用

rnorm

（或

runif

等）在 R 中生成随机数时，它们很少具有精确的均值和 SD 作为它们采样的分布。有没有简单的一两行代码可以帮我做到这一点？作为初步解决方案，我创建了这个函数，但它似乎应该是 R 或某个包的原生函数。

# Draw sample from normal distribution with guaranteed fixed mean and sd
rnorm_fixed = function(n, mu=0, sigma=1) {
  x = rnorm(n)  # from standard normal distribution
  x = sigma * x / sd(x)  # scale to desired SD
  x = x - mean(x) + mu  # center around desired mean
  return(x)
}

举例说明：

x = rnorm(n=20, mean=5, sd=10)
mean(x)  # is e.g. 6.813...
sd(x)  # is e.g. 10.222...

x = rnorm_fixed(n=20, mean=5, sd=10)
mean(x)  # is 5
sd(x)  # is 10

我想要这个的原因是我在将其应用于实际数据之前调整对模拟数据的分析。这很好，因为通过模拟数据，我知道确切的属性（均值、标准差等），并且我避免了 p 值膨胀，因为我正在做推论统计。我问是否存在像这样简单的东西

rnorm(n=20, mean=5, sd=10, fixed=TRUE)

Answer 1

既然你要求一句一句：

rnorm2 <- function(n,mean,sd) { mean+sd*scale(rnorm(n)) }
r <- rnorm2(100,4,1)
mean(r)  ## 4
sd(r)    ## 1

Answer 2

MASS 包中的 mvrnorm() 函数可以做到这一点。

library(MASS)
#empirical=T forces mean and sd to be exact
x <- mvrnorm(n=20, mu=5, Sigma=10^2, empirical=T)
mean(x)
sd(x)
#empirical=F does not impose this constraint
x <- mvrnorm(n=20, mu=5, Sigma=10^2, empirical=F
mean(x)
sd(x)

Answer 3

这是对先前答案中建议的功能的改进，以便它符合OP具有“固定”参数的需要。

仍然排成一行;-)

rnorm. <- function(n=10, mean=0, sd=1, fixed=TRUE) { switch(fixed+1, rnorm(n, mean, sd), as.numeric(mean+sd*scale(rnorm(n)))) }
rnorm.() %>% {c(mean(.), sd(.))}
#### [1] 0 1
rnorm.(,,,F) %>% {c(mean(.), sd(.))}
#### [1] 0.1871827 0.8124567

我选择为每个参数输入默认值，并添加

as.numeric

步骤来删除

scale

函数生成的属性。

Answer 4

作为附加答案，应该注意的是，问题中有一点含糊之处。

样本具有特定样本均值和样本方差的条件可以通过不同的方式执行。简单地标准化结果是一种方法，但在这种情况下，我们不总是从原始分布中进行条件采样（对于正态分布来说这并不重要）。

[0,1]范围内均匀分布的样本可以很好地说明这一点。如果我们对 n 个点进行采样并移动和重新调整它们的均值 0.5 和方差 1/12，那么某些点可能会落在 [0,1] 范围之外。

例如，像 {0,0.01,0.02,0.03,1} 这样的样本将重新调整为 {0.361 0.368 0.374 0.381 1.016}。

几何视点

从几何角度来看，我们可以将具有固定均值和方差的大小为 n 的样本视为固定在 n-2 维球体上的 n 维空间中的点。

在图中，我们展示了大小为 n=3 的均匀样本。每个样本都可以在欧几里得空间中描述，其中每个体积元素都具有相同的密度。下面我们通过画几个例子来说明这一点。

在第一个散点图中，我们绘制了所有 2000 个样本，并用红色表示样本均值和方差接近的情况。

在第二个散点图中，我们只绘制了样本均值和方差接近 1/2 和 1/12 的情况，它们是 [0,1] 范围内均匀分布的均值和方差。

我希望这张图片能够清楚地表明，具有特定方差和均值的样本位于垂直于对角线的圆形子空间上。

因此，在均值和方差固定的约束下从分布中进行采样的一般方法是使用拒绝采样，其中提案是从该圆中抽取的。

所以这可能是这样的：

runif2 = function(n,a,b) {

  ### compute variance and mean of uniform distribution 
  mean = (a+b)/2
  variance = (b-a)^2 / 12

  ### rejection sampling 
  reject = TRUE
  while( reject ) {
    ### simulate a sample on the circle 
    x = MASS::mvrnorm(n, mean, variance, empirical = TRUE) 

    ## compute density in proposed point 
    density = prod(dunif(x,a,b))
    density_max = 1/abs(b-a)^n

    ## rejection decision based on likelihood/density
    z = runif(1,0,density_max)
    reject = (z > density) 
  }
  return(x)
}

上面的代码中使用了拒绝采样，但对于均匀分布，也可以检查样本是否超出范围，而不是显式计算概率密度。