I am not sure off hand what may be wrong with your code, but optimization is only locally convergent and the radius of convergence can be quite small indeed, with cryptic error messages.
-nbreg- gets starting values from an ordinary Poisson regression, so maybe doing a bivariate Poisson with your variables is the way to go. Then you could hold those values constant and grid search over the negbin dispersion parameters to get better starts for those values, just evaluating the likelihood function for fixed values of the dispersion parameters. Once you've determined a good neighborhood for them, you can do a simultaneous refinement of all parameters.
I would generate simulated data that fits your model and see how it goes from different starts around the known parameter values. Jumping straight to real data with all its inherent messiness strikes me as a bridge too far.
In sum, my intuition is that adding the dispersion parameters changed a relatively benign log-likelihood into something quite a bit more complex.
-----Original Message-----
From: "Tirthankar Chakravarty" <[email protected]>
To: [email protected]
Sent: 3/12/2009 9:47 PM
Subject: st: Bivariate Negative Binomial using -ml-
Hello,
I am trying to implement a bivariate negative binomial regression
model as reported on pg. 212 of Rainer Winkelmann's book "Econometric
Analysis of Count Data" (Springer 5th edition 2008).
For those who do not have access to the book, it is the same
likelihood per observation reported in his paper ("Seemingly Unrelated
Negative Binomial Regression", Oxford Bulletin of Economics and
Statistics, Vol. 62, Issue 4, pg. 553-560, 2000) in equation (2)
substituting the NB-I density for the Poisson (also see equation (6)).
I have been trying to implement it using -ml- and my code is as below.
The problem is that -ml search- fails to find feasible initial values.
As -ml search- suggests, this may be because my code is "broken". I
have been looking at it but cannot spot the problem. Any help
debugging will be greatly appreciated.
Thanks and best regards,
Tirthankar
####################### BEGIN CODE################################
*********************
** Multivariate Negative Binomial (Ch. 7 pg.212, Winkelmann 2008) using -ml-
*********************
cap program drop binbreg2
program define binbreg2
version 9.1
args lnf xb1 xb2 sigma loggamma
tempvar z1 z2 mu1 mu2 fnb1 fnb2 fnb suma
qui gen double `mu1' = exp(`xb1')
qui gen double `mu2' = exp(`xb2')
local s = min($ML_y1, $ML_y2)
qui gen double `suma' = 0
qui gen float `z1' = .
qui gen float `z2' = .
qui gen float `fnb1' = .
qui gen float `fnb2' = .
qui gen float `fnb' = .
forvalues i = 0/`s' {
qui replace `z1' = $ML_y1 - `i'
qui replace `z2' = $ML_y2 - `i'
qui replace `fnb1' = lngamma(`mu1'/`sigma' +
`z1')-lngamma(`mu1'/`sigma') - lngamma(`z1' + 1) ///
-(`mu1'/`sigma')*log(1 + `sigma') +
`z1'*log(`sigma'/(1+`sigma'))
qui replace `fnb2' = lngamma(`mu2'/`sigma' +
`z2')-lngamma(`mu2'/`sigma') - lngamma(`z2' + 1) ///
-(`mu2'/`sigma')*log(1 + `sigma') +
`z2'*log(`sigma'/(1+`sigma'))
qui replace `fnb' = lngamma(`loggamma'/`sigma' +
`i')-lngamma(`loggamma'/`sigma') - lngamma(`i' + 1) ///
-(`loggamma'/`sigma')*log(1 + `sigma') +
`i'*log(`sigma'/(1+`sigma'))
qui replace `suma' = `suma' + `fnb1' + `fnb2' + `fnb'
}
qui replace `lnf' = `suma' if `suma' !=0 // Tried different
versions of this.
// qui replace `lnf' = . if `suma' ==0
end
ml model lf binbreg2 (first: breadth = sizeband lnemp lnturnover
geo_national geo_europe geo_other prodinnov dturnover_std
dturnover_std_sq) ///
(second: depth= sizeband lnemp lnturnover geo_national geo_europe
geo_other prodinnov dturnover_std dturnover_std_sq) /lna /loggamma
ml check
ml max
####################### END CODE################################
--
To every ω-consistent recursive class κ of formulae there correspond
recursive class signs r, such that neither v Gen r nor Neg(v Gen r)
belongs to Flg(κ) (where v is the free variable of r).
*
* For searches and help try:
* http://www.stata.com/help.cgi?search
* http://www.stata.com/support/statalist/faq
* http://www.ats.ucla.edu/stat/stata/
*
* For searches and help try:
* http://www.stata.com/help.cgi?search
* http://www.stata.com/support/statalist/faq
* http://www.ats.ucla.edu/stat/stata/
