Simulação Bivariate Poisson

Exemplo cálculo de probabilidade com a Poisson bivariada

\(P(Y1 = y1, Y2 = y2) = e^{-(\lambda x_{1}+\lambda x_{2}+\lambda x_{3})} \displaystyle\frac{\lambda x_{1}^{y1}}{y1!}\frac{\lambda x_{2}^{y2}}{y2!} \sum_{i=0}^{min(y1,y2)}\binom{y1}{i}\binom{y2}{i}i!\left(\frac{\lambda x_{3}}{\lambda x_{1}\lambda x_{2}}\right)^{i}\)

X1, X2 and X3 are independent Poisson variables with parameters: \(\lambda x_{1}, \lambda x_{2}, \lambda x_{3}\)

E(Y1)=\(\lambda x_{1}+\lambda x_{3}\)

E(Y2)=\(\lambda x_{2}+\lambda x_{3}\)

covariância(Y1,Y2) = \(\lambda x_{3}\)

correlação (Y1,Y2) = \(\lambda x_{3} / \sqrt{\lambda x_{1}+\lambda x_{3})(\lambda x_{2}+\lambda x_{3})}\)

library(extraDistr)
dbvpois(7, 7, 1, 1, 1) # Probability mass function Pr(y1=7,y2=7, lamda_x1=1,lamda_x2=1,lamda_x3=1 )

[1] 0.0002566068

dbvpois(1:7, 1:7, 1, 1, 1) # Probability mass function Pr(y1=1 a 7,y2= 1 a 7, lamda_x1=1,lamda_x2=1,lamda_x3=1 )

[1] 0.0995741367 0.0871273696 0.0470211201 0.0180650995 0.0053451950 0.0012799233 0.0002566068

Simulação de valores da distribuição de Poisson bivariada

Vamos observar na simulação que a covariância(Y1,Y2) = \(\lambda x_{3}\)

Vamos observar na simulação que a correlação (Y1,Y2) = \(\lambda x_{3} / \sqrt{\lambda x_{1}+\lambda x_{3})(\lambda x_{2}+\lambda x_{3})}\)

Para comparar o GC 2C com o GC ACS seria adequado usar o exemplo 1 para o GC 2C e o exemplo 2 para o GC ACS?

Exemplo 1

\(\lambda x_{1}=1, \lambda x_{2}=1, \lambda x_{3}=1\)

x <- rbvpois(100000, 1, 1, 1) # Simula poisson bivariada  -> lamda_x1=1,lamda_x2=1,lamda_x3=1
table(x[,1], x[,2]) # tabula os resultados

    
         0     1     2     3     4     5     6     7     8     9    10    11
  0   5034  4906  2539   848   201    40     5     0     0     0     0     0
  1   5038 10021  7591  3313  1000   254    44     8     0     0     0     0
  2   2415  7478  8681  5272  2079   642   162    35     7     0     0     0
  3    775  3302  5358  4766  2512   892   241    62    11     3     0     1
  4    209  1024  2212  2527  1839   903   309    98    10     1     0     0
  5     48   249   656   942   839   570   218    77    19     5     2     0
  6     14    49   142   267   304   224   135    52    28     5     0     0
  7      3     9    28    48   113    80    62    24    12     3     0     0
  8      1     0     3    14    12    24    21    10     3     2     0     0
  9      0     0     0     2     1     2     1     1     3     1     0     0
  10     0     0     0     1     0     0     1     1     1     0     0     0

round(prop.table(table(x[,1], x[,2])),3) # tabula a proporção dos resultados

    
         0     1     2     3     4     5     6     7     8     9    10    11
  0  0.050 0.049 0.025 0.008 0.002 0.000 0.000 0.000 0.000 0.000 0.000 0.000
  1  0.050 0.100 0.076 0.033 0.010 0.003 0.000 0.000 0.000 0.000 0.000 0.000
  2  0.024 0.075 0.087 0.053 0.021 0.006 0.002 0.000 0.000 0.000 0.000 0.000
  3  0.008 0.033 0.054 0.048 0.025 0.009 0.002 0.001 0.000 0.000 0.000 0.000
  4  0.002 0.010 0.022 0.025 0.018 0.009 0.003 0.001 0.000 0.000 0.000 0.000
  5  0.000 0.002 0.007 0.009 0.008 0.006 0.002 0.001 0.000 0.000 0.000 0.000
  6  0.000 0.000 0.001 0.003 0.003 0.002 0.001 0.001 0.000 0.000 0.000 0.000
  7  0.000 0.000 0.000 0.000 0.001 0.001 0.001 0.000 0.000 0.000 0.000 0.000
  8  0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
  9  0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
  10 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000

colMeans(x) # médias

[1] 2.00086 1.99946

cov(x[,1], x[,2]) # covariância = lamda_x3=1 

[1] 1.014821

cor(x[,1], x[,2]) #correlação = 1/raiz((1+1)x(1+1))

[1] 0.5046778

image(prop.table(table(x[,1], x[,2])))

Exemplo 2

\(\lambda x_{1}=2, \lambda x_{2}=2, \lambda x_{3}=0\)

y <- rbvpois(100000, 2, 2, 0) # Simula poisson bivariada  -> lamda_x1=2,lamda_x2=2,lamda_x3=0
table(y[,1], y[,2]) # tabula os resultados

    
        0    1    2    3    4    5    6    7    8    9   10
  0  1868 3690 3702 2399 1194  510  165   52   11    4    0
  1  3672 7333 7349 4947 2375  974  283   83   26    9    1
  2  3660 7280 7460 4926 2414  953  309   91   30    8    1
  3  2397 4767 4829 3277 1603  688  221   73   17    3    1
  4  1201 2445 2514 1595  868  322  108   32    8    3    1
  5   478  980  957  657  354  133   52   13    4    1    0
  6   164  310  312  199  103   36   17    6    0    0    0
  7    42  103  101   57   39   12    4    1    1    0    0
  8    13   25   18   17    7    1    0    1    0    0    0
  9     3    6    9    5    3    2    0    0    0    0    0
  10    0    2    0    0    0    0    0    0    0    0    0

round(prop.table(table(y[,1], y[,2])),3) # tabula a proporção dos resultados

    
         0     1     2     3     4     5     6     7     8     9    10
  0  0.019 0.037 0.037 0.024 0.012 0.005 0.002 0.001 0.000 0.000 0.000
  1  0.037 0.073 0.073 0.049 0.024 0.010 0.003 0.001 0.000 0.000 0.000
  2  0.037 0.073 0.075 0.049 0.024 0.010 0.003 0.001 0.000 0.000 0.000
  3  0.024 0.048 0.048 0.033 0.016 0.007 0.002 0.001 0.000 0.000 0.000
  4  0.012 0.024 0.025 0.016 0.009 0.003 0.001 0.000 0.000 0.000 0.000
  5  0.005 0.010 0.010 0.007 0.004 0.001 0.001 0.000 0.000 0.000 0.000
  6  0.002 0.003 0.003 0.002 0.001 0.000 0.000 0.000 0.000 0.000 0.000
  7  0.000 0.001 0.001 0.001 0.000 0.000 0.000 0.000 0.000 0.000 0.000
  8  0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
  9  0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
  10 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000

colMeans(x) # médias

[1] 2.00086 1.99946

cov(y[,1], y[,2]) # covariância = lamda_x3=1 

[1] 0.01368324

cor(y[,1], y[,2]) #correlação = 1/raiz((1+1)x(1+1))

[1] 0.006838898

image(prop.table(table(y[,1], y[,2])))