[Date Prev][Date Next][Thread Prev][Thread Next][Date index][Thread index]

From |
"David M. Drukker, Stata Corp" <ddrukker@stata.com> |

To |
statalist@hsphsun2.harvard.edu |

Subject |
Re: st: Random effects probit |

Date |
Wed, 18 Sep 2002 18:57:28 -0500 |

In the on going discussion about -xtprobit- Wiji Arulampalam <Wiji.Arulampalam@warwick.ac.uk> has raised several interesting questions. Let's begin with his comment that > I am dealing with models which have many variables . The models are > estimated using very large datasets. This is the first time I had > used xtprobit with lagged dependent variable (LDV) as a regressor. > When I didn't get any std. errors for rho or ln(sig2) but got proper > std. errors for the other coefficient estimates and when these > coefficient estimates were the same as the ordinary probit, I > automatically assumed that there was no unobserved heterogeneity! > Note, this problem only occurs when I have an L DV. But my > conclusions were obviously wrong - this I found out when I checked > it on Limdep. In the past I have always used Limdep for RE probit > and I have had no problems at all. The only reason I used xtprobit > was because of the size of the data set and model and I thought that > it would be easier given that my data was stored as a stata file. While there is some intuition here, the conclusion that there is no individual level heterogeneity does not hold. The reason is that you cannot compare two sets of parameter estimates when one of them does not come from a converged optimization procedure. Missing standard errors indicate that the esimated variance-covariance (VCE) matrix is not of full rank. After an ML procedure without constraints, an estimated VCE of less than full rank indicates that the optimization has not converged. So those unconverged -xtprobit- estimates with missing standard errors cannot be compared to the -probit- estimates. > It seems a bit too much to ask us to try different starting values > to check for sensitivity of the results when one is estimating a > very large model using a lot of data. How is one supposed to decide > whether the fact that the program gave the same answer as the > ordinary probit was because of unstable optimisation or whether it > was because of no unobserved heterogeneity? How would I > discriminate between these two? There are two points here. One is that is that the procedure Weihau and I recomended yesterday is going to be quite time consuming with large datasets. This is true. However, the problem that Wiji is having with missing standard errors is due to the fact that the current combination of optimization method and aproximation to the integral method used in -xtprobit- is not sufficiently robust. While the method of using a grid of starting values will be slow, it will probably overcome this lack-or-robustness for most datasets. Second, Wiji has asked how can one distinguish between unconverged results and no unobserved heterogeneity. In fact, he refined this question into two sub-questions in a subsequent post. > 1. If stata does not move from the rho=0.0 case and produces missing > std errors for rho and lnsig2, how do I know whether this is an > unstable numerical problem or that rho is genuinely zero? > 2. How come stata gives proper std errors for the other coefficient > estimates when it is also giving missing std errors for rho and > lnsig2? Allow to me provide a short and long version of my answer. Short version: If there are missing standard errors for any of the estimated coefficients from -xtprobit- and no constraints have been applied, treat the results as if they were from an a non-coverged maximization procedure. Try specifying different starting values. I continue to recommend the procedure that Weihau and I described yesterday. Second, if there are no missing standard errors and the results pass -quadchk-, then use a likelihood ratio test for inference with respect to the null hypothesis that there is no unobserved individual level heterogeneity. When no starting values are specified this test is reported at the bottom of the output. When starting values are specified, it must be calculated by hand. Below in my long version I give an example of how to perform this calculation. Long version: In this version, let's partition Wiji's questions into four parts. i) How can I distinguish converged from unconverged results? ii) If I obtain unconverged results with the default starting values, how can I try and get converged values. iii) If I have converged results from -xtprobit-, how can I tell if they are stable. iv) If I have converged, stable, results, how can test the null hypothesis of no unobserved heterogeneity. Let's begin with question i) and let's use simulated data. In the next several lines of output, I generate data from a simple cross-sectional -probit- model, although I do not provide an -id- variable for the panel routine to use. . set obs 500 obs was 0, now 500 . gen id = _n . expand 10 (4500 observations created) . set seed 123456 . gen double x1 = invchi2(2,uniform()) - 2 . gen double x2 = invchi2(2,uniform()) - 2 . gen double e = invnorm(uniform()) . gen double ystar = 3*(1 + x1 + x2)+ e . gen y = (ystar>0) Now, let's use -probit- to estimate the coefficients of this model. . probit y x1 x2 Iteration 0: log likelihood = -3439.6806 Iteration 1: log likelihood = -1820.8642 Iteration 2: log likelihood = -1136.0018 Iteration 3: log likelihood = -768.83952 Iteration 4: log likelihood = -588.73511 Iteration 5: log likelihood = -511.98949 Iteration 6: log likelihood = -491.23031 Iteration 7: log likelihood = -489.35275 Iteration 8: log likelihood = -489.33515 Probit estimates Number of obs = 5000 LR chi2(2) = 5900.69 Prob > chi2 = 0.0000 Log likelihood = -489.33515 Pseudo R2 = 0.8577 ------------------------------------------------------------------------------ y | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- x1 | 3.040576 .1270333 23.94 0.000 2.791595 3.289556 x2 | 2.991092 .1270602 23.54 0.000 2.742059 3.240125 _cons | 2.970903 .1338831 22.19 0.000 2.708496 3.233309 ------------------------------------------------------------------------------ note: 469 failures and 1394 successes completely determined. As one would expect, the -probit- estimates are very close to their true values of 3. Since we will need it below, let's save of the value of the log-likelihood. . scalar ll0 = e(ll) We will also find these estimated coefficients useful, so let's save them off as well. . mat b =e(b) Now, let's try to estimate the paramaters of an -xtprobit- model on this data. . xtprobit y x1 x2 , i(id) Fitting comparison model: Iteration 0: log likelihood = -3439.6806 Iteration 1: log likelihood = -1820.8642 Iteration 2: log likelihood = -1136.0018 Iteration 3: log likelihood = -768.83952 Iteration 4: log likelihood = -588.73511 Iteration 5: log likelihood = -511.98949 Iteration 6: log likelihood = -491.23031 Iteration 7: log likelihood = -489.35275 Iteration 8: log likelihood = -489.33515 Fitting full model: rho = 0.0 log likelihood = -489.33515 rho = 0.1 log likelihood = -490.43994 Iteration 0: log likelihood = -489.33514 Random-effects probit Number of obs = 5000 Group variable (i) : id Number of groups = 500 Random effects u_i ~ Gaussian Obs per group: min = 10 avg = 10.0 max = 10 Wald chi2(2) = 589.48 Log likelihood = -489.33514 Prob > chi2 = 0.0000 ------------------------------------------------------------------------------ y | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- x1 | 3.040577 .1270448 23.93 0.000 2.791574 3.28958 x2 | 2.991093 .1270721 23.54 0.000 2.742036 3.24015 _cons | 2.970904 .1338946 22.19 0.000 2.708475 3.233332 -------------+---------------------------------------------------------------- /lnsig2u | -14 . . . -------------+---------------------------------------------------------------- sigma_u | .0009119 . . . rho | 8.32e-07 . . . ------------------------------------------------------------------------------ Likelihood ratio test of rho=0: chibar2(01) = 0.00 Prob >= chibar2 = 1.000 Using the default starting values, we obtain results with missing standard errors. These are not converged results. So I make no attempt to interpret them. At this point, we move on to the second question: ii) If I obtain unconverged results with the default starting values, how can I try and get converged values. The answer is to try different starting values. I recomend using the -probit- estimates for the coefficients and a grid of points for rho. In my search for converged results I am going to specify some different starting values. And then use -xtprobit- in another attempt to obtain converged results. . mat b0 = b, (ln(.1/(1-.1))) . xtprobit y x1 x2 , i(id) from(b0, copy) Iteration 0: log likelihood = -490.43993 Iteration 1: log likelihood = -489.40071 Iteration 2: log likelihood = -489.25836 Iteration 3: log likelihood = -489.22553 Iteration 4: log likelihood = -489.22079 Iteration 5: log likelihood = -489.22057 Iteration 6: log likelihood = -489.22057 Random-effects probit Number of obs = 5000 Group variable (i) : id Number of groups = 500 Random effects u_i ~ Gaussian Obs per group: min = 10 avg = 10.0 max = 10 Wald chi2(2) = 363.32 Log likelihood = -489.22057 Prob > chi2 = 0.0000 ------------------------------------------------------------------------------ y | Coef. Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- x1 | 3.088922 .16568 18.64 0.000 2.764195 3.413649 x2 | 3.035387 .1604369 18.92 0.000 2.720936 3.349837 _cons | 3.015358 .1663938 18.12 0.000 2.689232 3.341484 -------------+---------------------------------------------------------------- /lnsig2u | -3.420701 2.189146 -7.711349 .8699469 -------------+---------------------------------------------------------------- sigma_u | .1808024 .1979015 .0211593 1.544922 rho | .0316547 .0671033 .0004475 .7047346 ------------------------------------------------------------------------------ Since there are no missing values for the standard errors, I have converged results. Unfortunately, specifying starting values prevented the likelihood ratio test from being automatically computed. Since I will have to do it by hand, I save off the value of the log-likelihood. . scalar ll1 = e(ll) Before I proceed with interpreting the results, I need to check that the converged results I have found do not vary unacceptably with the number of quadrature points used in aproximating the integral. In other words, I need to answer question iii). The answer is that I use -quadchk-. -quadchk- will reestimate the model for me with plus and minus 4 the number of quadrature points that I currently have. (You can change these numbers, see [R] quadchk.) If the estimates do not vary by much, then they are stable. . quadchk, noout Refitting model quad() = 8 Refitting model quad() = 16 Quadrature check Fitted Comparison Comparison quadrature quadrature quadrature 12 points 8 points 16 points ----------------------------------------------------- Log -489.22057 -489.22057 -489.22057 likelihood 3.968e-11 1.717e-11 Difference -8.110e-14 -3.509e-14 Relative difference ----------------------------------------------------- y: 3.0889218 3.0889218 3.0889218 x1 2.911e-11 -8.053e-12 Difference 9.424e-12 -2.607e-12 Relative difference ----------------------------------------------------- y: 3.0353869 3.0353869 3.0353869 x2 2.711e-11 -7.379e-12 Difference 8.931e-12 -2.431e-12 Relative difference ----------------------------------------------------- y: 3.0153583 3.0153583 3.0153583 _cons 2.665e-11 -7.402e-12 Difference 8.838e-12 -2.455e-12 Relative difference ----------------------------------------------------- lnsig2u: -3.420701 -3.420701 -3.420701 _cons 5.401e-10 -2.439e-10 Difference -1.579e-10 7.132e-11 Relative difference ----------------------------------------------------- Since the results are nearly invariant to the number of quadrature points, I conclude that they are stable and I can confidently proceed with interpretation. Now, let's look into quesiton iv). The test of the null that there is no heterogeneity is a test of the null hypothesis that sigma_u = 0. This is a test on the boundary of the parameter space for sigma_u. (See help j_chibar for more on this topic and a reference.) The computation of the test statisitic remains the same, but its asymptotic distribution is different than that of a stardard likelihood ratio test. Instead of converging to chi-squared with one degree of freedom it converges to another distribution which is basically .5*(chi-squared with one degree of freedom). Let's compute the value of the likelihood ratio test and its p-value. . scalar lr = 2*(ll1-ll0) . scalar p = .5*chi2tail(1,lr) And display the results: . di "hand lr is " lr " with p-value " p hand lr is .22917373 with p-value .31606859 This result indicates that we fail to reject the null hypothesis of no heterogeneity, or equivalently, that sigma_u =0 I hope that this helps. David ddrukker@stata.com * * For searches and help try: * http://www.stata.com/support/faqs/res/findit.html * http://www.stata.com/support/statalist/faq * http://www.ats.ucla.edu/stat/stata/

- Prev by Date:
**st: Surveys manipulation** - Next by Date:
**st: VARLAG help** - Previous by thread:
**Re: st: Random effects probit** - Next by thread:
**Re: st: Random effects probit** - Index(es):

© Copyright 1996–2014 StataCorp LP | Terms of use | Privacy | Contact us | What's new | Site index |