MajsoulPaipuAnalyzer 安定段位置题信区间计算

安定段位置题信区间计算问题 t-分布用于“用于根据小样本来估计呈正态分布且方差未知的总体的均值” 见百度百科：t分布 pt 在计算顺位后，我认为不是正态分布的（我猜测，因为排名有四个名次，而会有四个峰）因此置信区间不能准确反映安定段位的准确程度。

May 14 '19 08:05 kynowoo

你说的有道理，我模拟了数据以后确认了在t较小的时候很明显不是规则分布而是分布的叠加之类的形状。

但是形状也不是四个峰的，我个人猜想是123位贡献上差别不大只是加的PT有差别；四位次数会产生比较大的影响，因此并没有明显的四个峰。在N=10且安定段位可计算时目测是两个峰；而在N=100及以上时基本就看不出峰了，但是此时形状看起来有点像泊松分布；N=1000时就基本是正态分布了。

总之，在N较大时不论是安定段位对应的分布还是t分布都趋向于正态分布。在N较小时用t分布应该是不对的，但是我也不知道该用什么分布。。考虑到本来场次少时置信区间就没啥参考价值，以及和模拟结果比至少我目前结果的均值和方差都和真实分布的基本吻合，在知道真正分布前准备就先这么放着了。。

如果有安定段位遵循的分布的结论；或是更好的近似分布想法欢迎提供。

May 14 '19 10:05 zyr17

安定段位文档里提到的 $g(\hat{p}) = (A_1X_1+A_2X_2+A_3X_3)/X_4$ 是 $g(p) = (A_1p_1+A_2p_2+A_3p_3)/p_4$ 的MLE，但是是有偏估计，这是对的。安定段位的有偏性一节里，有小typo，在计算概率和期望的时候 $X_4$ 应该是 $p_4$。

但是目前的置信区间的计算方法是不准确的，文档里记录的现在的方法实际上在N小(<1000)和N大的时候(>10000)的时候都表现的不好。我尝试了两种理论上合理的计算置信区间的方法，用R生成并验证了他们在N不同的时候的10000次的模拟结果： N = 500，模拟：287.3585，446.6667 当前：303.3648，345.8318 方法1：293.0006，456.3935 方法2：281.8827，367.3140

N = 1000, 模拟：296.1000， 402.2785 当前：291.9306， 350.6110 方法1：298.6237，404.7964 方法2：291.6129，350.9287

N = 10000，模拟：310.2361，340.4440 当前：227.2820，408.6133 方法1：310.5515，340.6261 方法2：308.7354，327.1598

parametric bootstrap，直接用 $(X_1/N, X2/N, X3/N, X4/N)$ 作为多项分布的参数去生成 $m$ 个样本，计算这 $m$ 次的安定段位并且排序，取 $round(0.025m)$ 和 $round(0.975m)$ 处的值作为置信区间的端点（需要花一点点时间计算，但是很推荐）
delta method，核心是直接估计 $g(\hat{p})$ 的方差，再用正态分布去近似找到置信区间，也因为这个原因，其在N小的时候表现不好，计算式如下：

CI_delta <- function(X) {
  # X = (X1, X2, X3, X4) is the input vector
  # For simplicity, call MLE's p1,...,p4 in each run
  n_games <- sum(X)
  p1 <- X[1] / n_games
  p2 <- X[2] / n_games
  p3 <- X[3] / n_games
  p4 <- X[4] / n_games

  # Estimated gradient:
  est_gradient <- c(A1 / p4, A2 / p4, A3 / p4, -(A1 * p1 + A2 * p2 + A3 * p3) / p4^2)

  # Estimated var-cov matrix:
  est_var <- 1 / n_games^2 *
    rbind(
      c(n_games * p1 * (1 - p1), -n_games * p1 * p2, -n_games * p1 * p3, -n_games * p1 * p4),
      c(-n_games * p2 * p1, n_games * p2 * (1 - p2), -n_games * p2 * p3, -n_games * p2 * p4),
      c(-n_games * p3 * p1, -n_games * p3 * p2, n_games * p3 * (1 - p3), -n_games * p3 * p4),
      c(-n_games * p4 * p1, -n_games * p4 * p2, -n_games * p4 * p3, n_games * p4 * (1 - p4))
    )
  est_sd_theta <- sqrt(t(est_gradient) %*% est_var %*% est_gradient)
  return(c(
    estimator(X) - pnorm(0.975) * est_sd_theta,
    estimator(X) + pnorm(0.975) * est_sd_theta
  ))
}

Jan 14 '23 03:01 yuanhang0

关于把这个估计从有偏改成无偏的思路是计算 $E[g(\hat{p})] =E[ (A_1X_1+A_2X_2+A_3X_3)/X_4]$ 然后rescale到 $g(p) = (A_1p_1+A_2p_2+A_3p_3)/p_4$……不过计算有点繁琐，可能并不值得探索……有时间让wolfram alpha算一下好了。

Jan 14 '23 03:01 yuanhang0

@yuanhang0 感谢您的建议，确实直接采样是最方便的方法了，我自己拟合公式也是和采样结果拟合，当时可能考虑不想每次得出不一样的结果才找了个差不多的公式随便用一下。。请问你这里模拟时的次数和方法1中的m都取了多少？我之后考虑段位计算时固定随机种子用采样的方法实现。方法2的话，我不会rust所以看不太懂这段代码。。只不过用1就没什么问题了，我之后再慢慢研究吧，谢谢！

Jan 17 '23 15:01 zyr17

@yuanhang0 感谢您的建议，确实直接采样是最方便的方法了，我自己拟合公式也是和采样结果拟合，当时可能考虑不想每次得出不一样的结果才找了个差不多的公式随便用一下。。请问你这里模拟时的次数和方法1中的m都取了多少？我之后考虑段位计算时固定随机种子用采样的方法实现。方法2的话，我不会rust所以看不太懂这段代码。。只不过用1就没什么问题了，我之后再慢慢研究吧，谢谢！

不客气！模拟的次数我都用了10000次，方法1里的m我用了1000，因为我发现1000的时候表现就已经不错了而且基本不占计算时间。方法2是R的代码，我是用delta method去估计了这个东西的方差，样本小的时候表现不如1，在场次N=10000的时候才接近，实现1就可以了！

Jan 17 '23 21:01 yuanhang0

MajsoulPaipuAnalyzer MajsoulPaipuAnalyzer copied to clipboard

安定段位置题信区间计算

MajsoulPaipuAnalyzer
MajsoulPaipuAnalyzer copied to clipboard