r - 估计泊松分布

Question

我有一个图表，我计算了度数分布和度数，如下所示:

library(igraph) # for these two functions

dd <- degree_distribution(graph) 
d <- degree(graph)

由此，我估计到幂律，看看我的分布是否遵循“幂律”:

degree = 1:max(d)
probability = dd[-1]

nonzero.position = which(probability != 0)
probability = probability[nonzero.position]
degree = degree[nonzero.position]

reg = lm(log(probability) ~ log(degree))
cozf = coef(reg)

power.law.fit = function(x) exp(cozf[[1]] + cozf[[2]] * log(x))

从那以后，我使用 ggplot2 绘制了点和幂律。 .
结果如下图:

df <- data.frame(x = degree, y = probability)
  print(
      ggplot(df, aes(x,y,colour="Distribuição"))+
        geom_point(shape = 4) +
        stat_function(fun = power.law.fit, geom = "line", aes(colour="Power Law"))+

        labs(title = "Grafo", subtitle = "Distribuição dos Graus",
             x="K", y="P(k)", colour="Legenda")+
        scale_color_brewer(palette="Dark2")
  )

如您所见，我的分布不遵循幂律!我想估计泊松分布并在同一张图上绘图。
尽管我不确定我的分布不遵循(或遵循)泊松，但我想将它与幂律一起绘制。我不知道如何从数据中估计这个分布(泊松)，并计算平均度数。

谁能帮我？

The graph used to calculate the distribution and the degree is very large (700 thousand vertices), so I did not put the data of the graphs. The explanation of the answer can be based on any graph.

Answer 1

来自 ?dpois :

The Poisson distribution has density

p(x) = λ^x exp(-λ)/x!

for x = 0, 1, 2, … . The mean and variance are E(X) = Var(X) = λ.

因此，我将使用 secret lambda 生成一些虚拟数据:

mysecret <- ####

x <- data.frame(xes = rpois(50, mysecret))

> x$xes
 [1] 0 2 2 1 1 4 1 1 0 2 2 2 1 0 0 1 2 3 2 4 2 1 0 3 2 1 3 1 2 1 5 0 2 3 2 1 0 1 2 3 0 1 2 2 0 3 2 2 2 3


> mean(x$xes)
[1] 1.66
> var(x$xes)
[1] 1.371837

所以我的 secret lambda 的两个很好的猜测是 1.66 和 1.37。让我们试试看:

library(ggplot2)
ggplot(x, aes(xes)) + 
  geom_histogram(aes(y = ..density.., color = "Raw data"), 
                 fill = "white", binwidth = 1, center = 0, size = 1.5) +
  stat_summary(fun.y = dpois, aes(x = xes, y = xes, color = "Density based on E(X)"), 
               fun.args = list(lambda = 1.66), geom = "line", size = 1.5) +
  stat_summary(fun.y = dpois, aes(x = xes, y = xes, color = "Density based on Var(X)"), 
               fun.args = list(lambda = 1.37), geom = "line", size = 1.5)

他们俩都不错。你不能真正使用内置的 stat_function或 geom_density用于生成这些，因为泊松分布仅针对整数定义。不过，直方图和汇总函数运行良好，因为它们仅在数据点本身进行估计，而不是内插。

如果您想了解更多详细信息，可以使用 MASS包裹:

MASS::fitdistr(x$xes, dpois, start = list(lambda = 1))

    lambda  
  1.6601563 
 (0.1822258)

因此，让我们尝试从中构建:

library(dplyr)
df <- data_frame(xes = seq.int(max(x$xes)+1)-1,
                 dens.m = dpois(xes, 1.66),
                 dens.u = dpois(xes, 1.66+0.18),
                 dens.l = dpois(xes, 1.66-0.18))

> df
# A tibble: 6 x 4
    xes     dens.m     dens.u     dens.l
  <dbl>      <dbl>      <dbl>      <dbl>
1     0 0.19013898 0.15881743 0.22763769
2     1 0.31563071 0.29222406 0.33690378
3     2 0.26197349 0.26884614 0.24930880
4     3 0.14495866 0.16489230 0.12299234
5     4 0.06015785 0.07585046 0.04550717
6     5 0.01997240 0.02791297 0.01347012

ggplot(x, aes(xes)) + 
  geom_histogram(aes(y = ..density..), color = "black",
                 fill = "white", binwidth = 1, center = 0, size = 1.5) +
  geom_ribbon(data = df, aes(xes, ymin = dens.l, ymax = dens.u), fill = "grey50", alpha = 0.5) +
  geom_line(data = df, aes(xes, dens.m, color = "Based on E(X)\n+/-1 SD of lambda"), size = 1.5)

基于这两种方法和视觉解释，您应该可以轻松地说 λ = 1.66+/-0.18。

作为引用，我的 secret 初始值为 1.5。