Intermediate Econometrics

• Cases where the value of the response is a count, that is a non-negative integer
• Examples: annual doctor visits, number of cigarettes smoked daily, number of trips taken by members of households
  
• Drawbacks of fitting a linear model:
  • the integer nature of the response is not taken into account
  • the fitted model can predict negative values for some values of the covariates
  • if the distribution of the response is asymmetric, the logarithm transformation can’t be used in the common case where the response is 0 for a subset of observations

• Count data can be considered as coming from a random variable that follows a Poisson distribution, which has a unique parameter, this parameter being the mean and the variance of the series

\[Y_1, \ldots, Y_n \stackrel{i.i.d}{\sim} Poisson(\lambda), \quad \lambda > 0\] \[P(Y_i=y) = \frac{e^{-\lambda}\lambda^y}{y!}, \quad y\in \{0,1,2, \ldots\}\] \[\mathbf{E}(Y)=\lambda \quad \mathbf{V}(Y)=\lambda\]

• The ML estimator of this parameter is the sample mean \(\hat \lambda_MV = \bar Y\)

• Count data often take small values

• While comparing real count data to Poisson distribution, one often faces two problems:

  • overdispersion, which means that the variance is much greater than the mean
  • excess of zero, which means that the mean probability of a zero value computed using the Poisson distribution is often much less than the observed share of zero in the sample

• Overdispersion is important because we consider that the Poisson parameter is the same for every observation:

  • Adding covariates will give specific values of the Poisson parameter for every observations and will therefore reduce the (conditional) variance
    
  • Anyway, the overdispersion problem is often still present even after introducing covariates, because of unobserved heterogeneity

• The excess of zero problem is difficult to distinguish from overdispersion, as zero is an extreme value of the distribution and therefore, an excess of zero leads to overdispersion

  • The excess of zero can in part be explained by the fact that 0 is a special value that may not be correctly explained by the same process that explains strictly positive values

• Source: Terza (1998)
• Outcome: number of trips taken the day before the interview
• Let’s start by exploring the dataset and only considering the unconditional distribution (computing the first two moments and percentage of zero)

• This data exhibits important overdispersion, as the variance is about 5 times larger than the mean

• The predicted probability of 0 for the Poisson model is 0.01056, which is much lower than the actual percentage of 0 which is 18.5%

• There is therefore a large excess of 0, which also can be associated with overdispersion as 0 is a value far from the mean of 4.551

• Representation of the frequency of the data and probability mass of the fitted model

• For example, for \(Y=0\), the relative frequency is 0.18 and, as seen previously, the fitted probability is \(0.01056\) which implies

• The frequencies are clearly over-estimated for values of the response close to the mean (3, 4, 5, 6) and under-estimated for extreme values

• This is particularly the case for 0, which indicates that a specific problem of excess of zero may be present

• The model relies on the hypothesis that the response is a Poisson variable and therefore that the probability of a given value is

\[P(Y=y) = \frac{e^{-\mu}\mu^y}{y!}, \quad y\in \{0,1,2, \ldots\}\] having considered \(\mu = \lambda\)

• The Poisson distribution is defined by a unique parameter , which is the mean and the variance of the distribution. This is a striking feature of the Poisson distribution which differs for example from the normal distribution which has two separate position and dispersion parameters

• The Poisson probability can also be written as:

\[ \mbox{P}(Y = y, \mu) = e^{y \ln \mu-\mu- \ln y !} \]

• So, the Poisson model belongs to the generalized linear model family, with:

  • \(\theta = \ln \mu\)
  • \(b(\theta) = \mu = e^ \theta\)
  • \(\phi = 1\)
  • \(c(y, \phi) = - \ln y !\)

• The link defines the relation between the linear predictor \(\eta\) and the parameter of the distribution \(\mu\): \(\eta = g(\mu)\). The canonical link is, for the Poisson model, the log link: \(\eta = \ln \mu\)

• For a given set of covariates, the Poisson parameter for the \(i\)-th observation is \(\mu_i = g ^ {-1}(\eta_i)\), with \(\eta_i = x^\top_i \beta\) the linear predictor

\[ \mu_i = \mbox{E}(Y_i \mid x_i) = e^{x^\top_i \beta} \]

and the log-likelihood function is:

\[ \ln L (\beta) = \sum_{i=1} ^ n \ln \frac{e^{-e^{x^\top_i \beta}} e ^{y_ix^\top_i \beta}}{y_i !} = \sum_{i=1} ^ n -e^{x^\top_i \beta} + y_i x^\top_i \beta - \ln y_i ! \]

• We can implement the log-likelihood and optimize it by using optimizer, or implementing the Newton-Raphson algorithm requiring the first and second order derivative (the latter is useful for obtaining the standard errors of the parameter estimator)

• The covariates are:
  • the share of trips for work or school (workschl)
  • the number of individuals in the household (size)
  • the distance to the central business district (dist)
  • a dummy (smsa) for large urban area
  • the number of full-time workers in a household (fulltime)
    
  • the distance from home to the nearest transit node (distnode)
  • household income divided by the median income of the census tract (realinc)
    
  • a dummy if the survey period is Saturday or Sunday (weekend)
  • a dummy for owning at least one car (car)

• The Poisson model is fitted using the glm function and setting the family argument to poisson
• For the sake of comparison, we also estimate a linear model using lm:

• An increase of 1 of realinc means that household’s income increases by an amount equal to the local median income
• The linear models indicates that the number of trips then increases by \(0.133\) trips
• The Poisson model indicates:
  • that the relative increase of trips \(dy /y\) is equal to \(0.019\), i.e., 2%
  • At the sample mean, it corresponds to an increase of \(0.019 \times 4.55 = 0.086\) trips

• If the interview concerns a weekend day, the OLS coefficient is \(-0.48\), which indicates that the number of trips taken during the weekend are about one-half less, compared to a weekday

• For such a dummy variable, denoting \(\beta_w\) and \(\tilde x\) a given vector of covariates such that \(x_w = 0\), \(\tilde x^\top_i\beta\) is the linear predictor for a weekday and \({\tilde x}^\top_{i}\beta + \beta_w\) the linear predictor for a week-end day

• Therefore, the relative difference \(\tau\) of the number of trips between a weekend day and a weekday is

\[ \tau = \frac{P(Y=k| weekend,rest) - P(Y=k| weekday,rest)}{P(Y=k| weekday,rest)}\]

\[\tau = \frac{e^{\tilde x^\top_i\beta + \beta_w}- e^{\tilde x^\top_i \beta}}{e^{\tilde x^\top_i \beta}} = e^{\beta_w} - 1 = e^{-0.0738} - 1 = - 0.071 = -7.1\% \]

• The previous expression indicates also that \(\beta_w = \ln (1 + \tau) \approx \tau\). For small values of \(\beta_w\), the coefficient can be interpreted as the relative difference of the response

• In terms of absolute value, at the sample mean, the absolute difference is \(-0.071 \times 4.55 = -0.323\)

• The matrix of the second derivatives is:

\[ \frac{\partial \ln L}{\partial \beta\partial \beta ^ \top}=-\sum_{i=1}^n e^{x^\top_i \beta}x_ix^\top_i =- X^\top diag(e^{x^\top_1 \beta}, \ldots, e^{x^\top_n \beta}) X \]

• The information matrix is the opposite of the hessian (as it doesn’t depend on \(y\), it is equal to its expected value)

• The estimated variance of the estimator is therefore:

\[ \hat{\mbox{V}}(\hat{\beta}) = \left(\sum_{i=1}^n e^{x^\top_i \hat \beta}x_ix^\top_i\right) ^ {-1} = \left(\sum_{i=1} ^ n \hat{\mu}_i x_i x^\top_i\right) ^ {-1} \]

• This estimator is consistent only if the distribution of the response is Poisson, and therefore if the conditional variance equals the conditional mean

• A more general estimator is based on the sandwich formula: \[ \left(\sum_{i=1} ^ n \hat{\mu}_i x_i x^\top_i\right) ^ {-1} \left(\sum_{i=1} ^ n (y_i - \hat{\mu}_i)^2 x_i x^\top_i\right) \left(\sum_{i=1} ^ n \hat{\mu}_i x_i x^\top_i\right) ^ {-1} \]

• One alternative is to assume that the variance is proportional to the mean: \(\mbox{V}(Y|x) = \phi \mbox{E}(Y|x)\). \(\phi\) can then be consistently estimated by:

\[ \hat{\phi} = \frac{1}{n}\sum_{i=1} ^ n \frac{(y_i - \hat{\mu}_i) ^ 2}{\mu_i} \]

which leads to the third estimator of the covariance matrix:

\[ \hat{\mbox{V}}(\hat{\beta}) = \hat{\phi}\left(\sum_{i=1} ^ n \hat{\mu}_i x_i x^\top_i\right) ^ {-1} \]

• The “quasi-Poisson” model is obtained by using the log link and the equation above for the variance. It leads to the same estimator as the Poisson model and a variance that is greater (respectively lower) than the one of the Poisson model if there is overdispersion (respectively underdispersion)

• The quasi-poisson model can be fitted using glm by setting the family argument to quasipoisson
• The coefficient estimates are the same of those obtained via Poisson

• Compared to the Poisson model, the standard deviations of the quasi-Poisson model are inflated by a factor which is the square root of the estimate of the variance-mean ratio

• This factor can be estimated using

\[ \hat{\mbox{V}}(\hat{\beta}) = \hat{\phi}\left(\sum_{i=1} ^ n \hat{\mu}_i x_i x^\top_i\right) ^ {-1} \]

where \[ \hat{\phi} = \frac{1}{n}\sum_{i=1} ^ n \frac{(y_i - \hat{\mu}_i) ^ 2}{\mu_i} \] which makes use of the residuals of the Poisson regression:

• The response residual is \(y_i - \hat{\mu}_i\)

• The Pearson residual are obtained by dividing the response residual by the standard deviation, which is \(\sqrt{\hat{\mu}_n}\)

• Then, the estimation of the variance-mean ratio is simply the sum of square of the Pearson residuals divided by the sample size (the number of degrees of freedom can be used instead)

• The sandwich estimator can be constructed by extracting from the fitted Poisson model the model matrix (\(X\)), the fitted values (\(\hat{\mu}_i\)) and the response residuals (\(y_i - \hat{\mu}_i\))

\[ \left(\sum_{i=1} ^ n \hat{\mu}_i x_i x_i\right) ^ {-1} \left(\sum_{i=1} ^ n (y_i - \hat{\mu}_i)^2 x_i x_i^\top\right) \left(\sum_{i=1} ^ n \hat{\mu}_i x_i x_i\right) ^ {-1} \]

• A more general expression for the variance is:

\[ \sigma_i ^ 2 = \mu_i + \alpha \mu_i ^ p \]

• We will focus on \(p = 1\) which implies that the variance is proportional to the expectation, \(\sigma_i ^ 2 = (1 + \alpha) \mu_i\),

• Tests for overdispersion can be done using one of the three test principles: Wald, likelihood ratio and score tests:

  • The first two require the estimation of a more general model that doesn’t impose the equality between the conditional mean and the conditional variance
  • The score test requires only the estimation of the constrained (Poisson) model: we consider an auxiliary regression.^[Score tests that are not obtained from auxiliary regressions approach

• The hypothesis of equidispersion (i.e., \(\alpha = 0\)) can be tested by regressing \(\frac{ (y_n - \hat{\mu}_n) ^ 2 - y_n}{\hat{\mu}_n}\) on a constant and comparing the Student statistic to the relevant critical value

• In both regressions the p-value is very close to zero so that the equidispersion hypothesis is strongly rejected

• The overdispersion phenomenon can also be represented by plotting the square of the residuals against the fitted values
• The cloud of points forms a sort of ‘fan’ opening upward meaning that the variance is greater than mean

• We can handle with overdispersion using a different model: using the negative binomial distribution to describe the distribution of \(Y\) conditional to \(x\)

• Under this model \[\mathbb{E}(Y_i |x) = \mu_i\] \[\mathbb{V}(Y_i |x) = \mu_i(1+\mu_i/\delta)\]

• Model fitting can be carried out using the glm.nb function

• In the previous model output summary the value of \(\theta\) corresponds to \(1/\delta\) and provide support for having a overdispersion problem
• The standard errors are much higher compared to the Poisson model, but remember that, when overdispersion is present, the standard errors of the Poisson model are downward biased
• Finally, the AIC clearly shows that the binomial negative model outperform the poisson one

• One of the problem trips data set is that the Poisson model is unable to model correctly the probability of \(Y = 0\)

• Different phenomena can explain zero values:

  • some people may never (or almost never) take a trip out of the house because, for example, of a bad physical condition,
  • some people may take a trip infrequently (for example, twice a week) so that a 0 value may be reported if the survey concerns a Tuesday as the individual has taken a trip on Monday and Thursday

• A viable option to consider (there exists also other alternatives) is the Hurdle model

  • a binomial model that explains the fact that the outcome is zero or positive
  • a model that explains the level of the outcome when it is positive

• The Hurdle model is general, meaning that it can be used for several type of outcomes

• In our case we are modelling \(P(Y=0)\) and the \(P(Y=k|Y>0)\) using a Poisson or Negative Binomial

• We can fit these models using the countreg R package. However, it is on the CRAN and must be install via github

• It seems that the Hurdle model is responsible for an improved model fit and the preferred distribution choice is the Negative binomial

• Denoting \(W\) the matrix of instruments, the IV estimator is \[\hat{\beta} = (X^\top P_W X) ^ {-1} (X^\top P_W y)\]

• While the GMM estimator is
  \[\hat{\beta} = \left(X^\top W \hat{S}^{-1} W ^ \top X\right)^{-1}X^\top W \hat{S}^{-1} W ^ \top y\] with

  • \[\hat{S} = \sum_i \hat{\epsilon}_n ^ 2 w_i w^\top_i\]
  • \(\hat{\epsilon}_i\) being the residual of a consistent estimation

• Source: Mullahy (1997)

• Estimate a demand function for cigarettes which depends on the stock of smoking habits

• This variable is quite similar to a lagged dependent variable and is likely to be endogenous as the unobservable determinants of current smoking behavior should be correlated with the unobservable determinants of past smoking behavior

• The data set, called cigmales, contains observations of 6160 males in 1979 and 1980 from the smoking supplement to the 1979 National Health Interview Survey

• The response cigarettes is the number of cigarettes smoked daily

• The covariates are:
  • the habit “stock” habit
  • the current state-level average per-pack price of cigarettes price
  • a dummy indicating whether there is in the state of residence a restriction on smoking in restaurants restaurant
  • the age age
  • the number of years of schooling educ and their squares
  • the number of family members famsize
  • a dummy race which indicates whether the individual is white or not
• The external instruments are:
  • cubic terms in age and educ and their interaction
  • the one-year lagged price of a pack of cigarettes lagprice
  • the number of years the state’s restaurant smoking restrictions had been in place

• Let’s start to padd the dataset with additional variables (square and cubic terms and interaction age and education)

• The starting point is a basic count model, i.e., a Poisson model with a log-link and robust standard errors

• The IV and the GMM estimators are estimated using the micsr::expreg function

• Its main argument is a two-part formula, the first part indicating the covariates and the second part the instruments
  

• The method argument can be set to "iv" or "gmm" to estimate respectively the instrumental variable and the general method of moments estimators

• For the two flavors of instrumental variable estimators, the coefficient of habit is about half of the one of the Poisson model, indicating that this covariate is positively correlated with the unobserved determinants of smoking

• Used in IV and GMM estimation when there are more instruments than endogenous variables

• Null hypothesis: all instruments are valid (they affect the dependent variable only through the endogenous regressors and are uncorrelated with the structural errors).

• Test statistic is based on the sample moments of the instruments

• Low p-value then instruments fail and the model can be misspecifid or the instruments are invalid

• The Sargan test suggest exogeneity of the instruments hypothesis is not rejected