R help - Oct 2008 - Negative Binomial Predictions

Good Day All,

I have a negative binomial model which I have developed using the MASS 
library.  I now would like to develop some predictions from it. 

Running the predict.glm (stats library) using type="response" gives me
a
non-integer value which was rather puzzling.  I would like to confirm 
that this is actually the mean predicted value of the probability mass 
function as opposed to the most likely value.  I am afraid that reading 
the help file for predict.glm either does not state this or I don't 
understand what it states (which of course could always be the case)

Given that this is a negative binomial model, the mean is often times to 
the right of the most likely value, so I'd like to ask how one would go 
about predicting the most likely value. 

Thanks in advance !

-- 

-Don 

Don Catanzaro, PhD                  Landscape Ecologist
dgcatanzaro at gmail.com               16144 Sigmond Lane
479-751-3616                        Lowell, AR 72745

On Wed, 1 Oct 2008, Donald Catanzaro, PhD wrote:
> Good Day All,
>
> I have a negative binomial model which I have developed using the MASS 
> library.  I now would like to develop some predictions from it. 
> Running the predict.glm (stats library) using type="response"
gives me a
> non-integer value which was rather puzzling.  I would like to confirm that 
> this is actually the mean predicted value of the probability mass function
as
> opposed to the most likely value.  I am afraid that reading the help file
for
> predict.glm either does not state this or I don't understand what it
states
> (which of course could always be the case)
It gives the conditional mean, i.e., the inverse link applied to the 
linear predictor x'b.
> Given that this is a negative binomial model, the mean is often times to
the
> right of the most likely value, so I'd like to ask how one would go
about
> predicting the most likely value.
You can compute the expected probabilities for each count using the 
predprob() method in package "pscl", or by hand using the dnbinom() 
function:

## load data and fit model
library("pscl")
data("bioChemists")
fm <- glm.nb(art ~ ., data = bioChemists)

## compute expected probabilites for y = 0
## via predprob()
pp <- predprob(fm)
pp[,1]
## by hand
dnbinom(0, mu = fitted(fm), size = fm$theta)

## compute the count with the highest probability for each observation
apply(pp, 1, which.max) - 1
## or the median
apply(pp, 1, function(x) max(which(cumsum(x) < 0.5)))

hth,
Z
> Thanks in advance !
>
> -- 
>
> -Don 
> Don Catanzaro, PhD                  Landscape Ecologist
> dgcatanzaro at gmail.com               16144 Sigmond Lane
> 479-751-3616                        Lowell, AR 72745
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
>