在EDX r Stats类中,我们被要求查看'6'''6'的比例在100个模具卷中滚动。然后,我们被要求滚动100次10,000次掷骰,以查看100次卷中平均值的标准偏差。
100-DIE卷的结果是预期的;大约0.1703左右(1/6 = 0.1666667)
但是,当我加载replicate()以投掷100骰子的套件10,000次以查看10,000个均值时,结果不是我所期望的。我看不到z得分= 2的范围之外的任何值:
set.seed(1)
# get mean of 100 dice rolls
mean100dice <- function(){
n=100
x <- replicate(n, sample(1:6, n, replace=TRUE), simplify='vector')
mean(x==6)
}
mean100dice() #these come out as expected
means10k <- replicate(10000, mean100dice(),simplify='vector')
p = 1/6
z = (means10k - p) / sqrt(p*(1-p)/n)
mean(z > 2) ## I expect this to be > 0
range(means10k) ## sanity check
> mean(z > 2)
[1] 0
> range(means10k)
[1] 0.1522 0.1806
在猜测中,您在计算z
时设置n <- 100
而不是n <- 10000
。
提供明确的变量名称是一个好主意,因此您不会混合。例如,您需要区分n_dice_rolls
和n_replicates
。
顺便说一句,您计算100个骰子卷的平均值的代码不正确。
sample(1:6, n, replace=TRUE)
卷n
骰子;您也不需要致电replicate()
。我认为您想要这样的东西。
roll_nd6 <- function(n_dice) {
sample(1:6, n_dice, replace = TRUE)
}
get_fraction_of_sixes_from_rolling_nd6 <- function(n_dice) {
mean(roll_nd6(n_dice) == 6L)
}
monte_carlo_simulate_get_fraction_of_sixes <- function(n_replications, n_dice) {
replicate(
n_replications,
get_fraction_of_sixes_from_rolling_nd6(n_dice),
simplify = "vector"
)
}
calc_z_score <- function(actual_p, expected_p) {
(actual_p - expected_p) /
sqrt(expected_p * (1 - expected_p) / length(actual_p))
}
actual_fraction_of_sixes <- monte_carlo_simulate_get_fraction_of_sixes(10000, 100)
z_scores <- calc_z_score(actual_fraction_of_sixes, 1 / 6)
您在mean100dice
中有一个错误:您采样了100个骰子,然后重复100次,因此实际上不是100骰子的平均值,而是100*100*100 = 10,000骰子。当然,平均水平将更接近P。