R 的 poweRlaw 包中连续与离散对数正态分布的似然函数差异

Question

我正在尝试使用 Colin Gillespie 在 R 中的 poweRlaw 包将对数正态分布拟合到一些计数数据。我知道对数正态分布是连续的，而计数数据是离散的，但是，包包含类和对数正态分布的连续和离散版本的方法。

当我拟合 xmin（忽略计数值的阈值）、log mean 和 log sd 参数和 bootstrap 结果以获得 p 值时，我得到向量内存耗尽错误。我发现当包内部函数 sample_p_helper 试图从拟合分布生成随机数时会发生这种情况。拟合的 log mean 和 log sd 参数太低以至于拒绝抽样算法试图生成字面上数十亿的数字以获得高于 xmin 的任何值，因此内存问题。

输入：

library(poweRlaw)

counts <- c(54, 64, 126, 161, 162, 278, 281, 293, 296, 302, 322, 348, 418, 511, 696, 793, 1894)

dist <- dislnorm$new(counts)      # Create discrete lnorm object
dist$setXmin(estimate_xmin(dist)) # Get xmin and parameters

bs <- bootstrap_p(dist) # Run bootstrapping

错误信息：

Expected total run time for 100 sims, using 1 threads is 24.6 seconds.
Error in checkForRemoteErrors(val) : 
  one node produced an error: vector memory exhausted (limit reached?)

那么问题就变成了为什么首先要拟合如此低且拟合不佳的 log mean 和 log sd 参数值。

我注意到如果我拟合对数正态分布的连续版本，错误不会发生并且参数值似乎更合理（事实上，p 值表明数据符合对数正态分布）：

dist_cont <- conlnorm$new(counts)
dist_cont$setXmin(estimate_xmin(dist_cont))

bs <- bootstrap_p(dist_cont)

bs

查看包的源代码，我注意到离散对数正态分布与连续对数正态分布的似然函数不同。具体来说就是计算联合概率的部分。

The continuous version 看起来和我期望的一样：

########################################################
#Log-likelihood 
########################################################
conlnorm_tail_ll = function(x, pars, xmin) {
  if(is.vector(pars)) pars = t(as.matrix(pars))
  n = length(x)
  joint_prob = colSums(apply(pars, 1, 
                             function(i) dlnorm(x, i[1], i[2], log=TRUE)))

  prob_over = apply(pars, 1, function(i) 
    plnorm(xmin, i[1], i[2], log.p=TRUE, lower.tail=FALSE))
  joint_prob - n*prob_over

}

然而，在the discrete version中，联合概率的计算方式不同：

########################################################
#Log-likelihood 
########################################################
dis_lnorm_tail_ll = function(xv, xf, pars, xmin) {
  if(is.vector(pars)) pars = t(as.matrix(pars))
  n = sum(xf)
  p = function(par) {
    m_log = par[1]; sd_log = par[2]
    plnorm(xv-0.5, m_log, sd_log, lower.tail=FALSE) - 
      plnorm(xv+0.5, m_log, sd_log, lower.tail=FALSE)
  }
  if(length(xv) == 1L) {
    joint_prob = sum(xf * log(apply(pars, 1, p)))
  } else {
    joint_prob = colSums(xf * log(apply(pars, 1, p)))
  }
  prob_over = apply(pars, 1, function(i) 
    plnorm(xmin-0.5, i[1], i[2], 
           lower.tail = FALSE, log.p = TRUE))

  return(joint_prob - n*prob_over)
}

指数分布的离散和连续实现之间存在类似差异，但离散和连续幂律分布之间没有差异。在连续版本中，joint_prob 是通过相对简单地调用 dlnorm 计算的，但离散版本改为调用 plnorm。此外，他们调用 plnorm 两次，首先是观察数据值 -0.5，然后是观察值 +0.5，然后从后者中减去前者。

所以，最后，我的问题是：

为什么poweRlaw在对数正态分布的离散实现中会这样计算联合概率？我敢肯定它是以这种方式写的是有原因的，这只是我的数学无知，但我真的不明白。
使用 poweRlaw 的连续对数正态分布是否安全，即使我的数据是离散的，因为它似乎工作得很好？
关于在尝试拟合离散对数正态分布时我的数据可能出了什么问题的任何其他线索？我在想某个地方可能存在缩放问题，但我很难理解它。
我的小得可笑的数据集能发挥作用吗？我正在尝试将分布拟合到仅高于 xmin 的 8 个值，我知道这对于最大似然可靠来说太少了。

谢谢你们陪我度过这漫长的post。我知道这既是一个统计问题，也是一个编码问题。非常感谢任何朝着正确方向有用的推动！干杯。

Answer 1

dlnorm()给出概率密度值。记住密度积分为一但总和不为一。因此，为了计算出离散分布，我们取整数两边的值。它们也将是一个标准化常数。对于 CTN 案例，对数似然只是 dlnorm() 的乘积，这样更容易也更快。
"Safe" 是一个很难定义的词。对于此数据，CTN 和离散在视觉上给出了相同的拟合。但两者都不太合适。
离散分布的估计参数值在非常极端的尾部给出截断的对数正态分布。在该地区模拟数据具有挑战性
是的，您的数据有问题。但这也是模型不起作用时的问题；）

R 的 poweRlaw 包中连续与离散对数正态分布的似然函数差异

Difference in likelihood functions for continuous vs discrete lognormal distributions in R's poweRlaw package

r

distribution

power-law

data-fitting