R glm（）の係数からポアソン分布 "lambda"を得るには

Question

私のRスクリプトは以下のglm（）coeffsを生成します。 ポアソンのラムダとは何ですか？これは私がディストリビューションを作成するために使用したものなので、〜3.0にする必要があります。

Call: 
glm(formula = h_counts ~ ., family = poisson(link = log), data = pois_ideal_data) 

Deviance Residuals: 
    Min  1Q Median  3Q  Max 
-22.726 -12.726 -8.624 6.405 18.515 

Coefficients: 
      Estimate Std. Error z value Pr(>|z|)  
(Intercept) 8.222532 0.015100 544.53 <2e-16 *** 
h_mids  -0.363560 0.004393 -82.75 <2e-16 *** 
--- 
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

(Dispersion parameter for poisson family taken to be 1) 

Null deviance: 11451.0 on 10 degrees of freedom 
Residual deviance: 1975.5 on 9 degrees of freedom 
AIC: 2059 

Number of Fisher Scoring iterations: 5 


random_pois = rpois(10000,3) 
h=hist(random_pois, breaks = 10) 
mean(random_pois) #verifying that the mean is close to 3. 
h_mids = h$mids 
h_counts = h$counts 
pois_ideal_data <- data.frame(h_mids, h_counts) 
pois_ideal_model <- glm(h_counts ~ ., pois_ideal_data, family=poisson(link=log)) 
summary_ideal=summary(pois_ideal_model) 
summary_ideal 

Answer 1

ここで何をしていますか？分布に合うようにglmを使用しましたか？

まあ、そうすることは不可能ではありませんが、これを介して行われます：

set.seed(0) 
x <- rpois(10000,3) 
fit <- glm(x ~ 1, family = poisson()) 

はつまり、私たちは、切片のみの回帰モデルを使用してデータをフィット。

fit$fitted[1] 
# 3.005 

これは同じです：あなたは集計やビニングデータにフィットポアソンをやろうとしているように見えます

mean(x) 
# 3.005 

Answer 2

。それはglmのことではありません。私は、缶詰めのデータにディストリビューションをフィッティングするための缶詰の方法をすばやく見ましたが、見つけられませんでした。 bdaパッケージの以前のバージョンはこれを提供しているかもしれませんが、今はそうではないようです。

rootで、(# counts)*prob(count|lambda)を計算し、optim()を使用して最小化する負の対数尤度関数を設定する必要があります。 bbmleパッケージを用いて以下に与えられたソリューションは、もう少し複雑なアップフロントですが、あなたは簡単に計算信頼区間などのようなメリットを追加できます。..データを設定し

：

set.seed(101) 
random_pois <- rpois(10000,3) 
tt <- table(random_pois) 
dd <- data.frame(counts=unname(c(tt)), 
       val=as.numeric(names(tt))) 

ここで私はtableを使用していますむしろhistより離散データのヒストグラムがbbmleで動作するようにビニングポアソンデータの密度関数を（ '設定し

（あなたは左閉鎖対右に注意する必要があるため、整数カットポイントを持つことがしばしば物事が混乱します）うるさいですので、フォーム最初の引数はxとする必要があり、引数はlogでなければなりません）。

library(bbmle) 
m1 <- mle2(counts~dpoisbin(val,exp(loglambda)), 
     data=dd, 
     start=list(loglambda=0)) 
all.equal(unname(exp(coef(m1))),mean(random_pois),tol=1e-6) ## TRUE 
exp(confint(m1)) 
## 2.5 % 97.5 % 
## 2.972047 3.040009 

：

dpoisbin <- function(x,val,lambda,log=FALSE) { 
    probs <- dpois(val,lambda,log=TRUE) 
    r <- sum(x*probs) 
    if (log) r else exp(r) 
} 

フィットラムダは、（リンクは負のラムダ値から数値的な問題/警告を防ぐことができますログ）