Notice: On April 23, 2014, Statalist moved from an email list to a forum, based at statalist.org.
From | Jorge Eduardo Pérez Pérez <jorge_perez@brown.edu> |
To | "statalist@hsphsun2.harvard.edu" <statalist@hsphsun2.harvard.edu> |
Subject | Re: st: Error in Poi's QUAIDS user-command |
Date | Fri, 18 Apr 2014 15:32:27 -0400 |
Which version of -quaids- are you using? Mine is version 1.1.0 and the code reads: if (quadratics == "") { bofp = exp(lnp*beta') lnpindex = lnp*alpha' :+ st_numscalar("e(anot)") for(i=1; i<=rows(lnpindex); ++i) { lnpindex[i] = lnpindex[i] + 0.5*lnp[i,.]*gamma*lnp[i,.]' } } -------------------------------------------- Jorge Eduardo Pérez Pérez Graduate Student Department of Economics Brown University On Fri, Apr 18, 2014 at 2:43 PM, Robin Winkler <robinwinkler@gmail.com> wrote: > Hi all, > > I believe I have found a serious error in the way in which the QUAIDS > user command written by Brian Poi computes elasticities. > > In the Mata code underlying the ado file - the relevant bit is pasted > below - the price index *lnpindex* used to compute the elasticities > does not include alpha 0, which is entered into the command suite as > *anot*. Hence, I have had the problem that my price index was > underestimated, and the magnitude of the computed expenditure > elasticities tend to be way overstated. > > Should Mr Poi read this, I would be grateful for a comment. > > Robin > > > > *! version 1.0.0 29dec2011 > /* > _quaids__utils.mata > Mata routines called by > nlsur__quaids.ado > quaids_estat.ado > quaids_p.ado > */ > > mata > > mata set matastrict on > > void _quaids__expelas(string scalar touses, > string scalar quadratics, > string scalar atmeans, > string scalar lnps, > string scalar lnexps, > real scalar ndemo, > string scalar demos, > string scalar expelass) > { > real scalar i > real vector alpha, beta, lambda, rho > real vector bofp, cofp, lnpindex, lnexp, mbar > real matrix gamma, eta > real matrix shares, expelas, lnp, demo > alpha = st_matrix("e(alpha)") > beta = st_matrix("e(beta)") > gamma = st_matrix("e(gamma)") > if (quadratics == "") { > lambda = st_matrix("e(lambda)") > } > else { > lambda = J(1, cols(beta), 0) > } > if (ndemo > 0) { > eta = st_matrix("e(eta)") > rho = st_matrix("e(rho)") > if (atmeans == "") { > st_view(demo=., ., demos, touses) > } > else { > demo = mean(st_data(., demos, touses)) > } > } > > if (atmeans == "") { > st_view(shares=., ., st_global("e(lhs)"), touses) > st_view(expelas=., ., expelass, touses) > st_view(lnp=., ., lnps, touses) > st_view(lnexp=., ., lnexps, touses) > } > else { > shares = mean(st_data(., st_global("e(lhs)"), touses)) > lnp = mean(st_data(., lnps, touses)) > lnexp = mean(st_data(., lnexps, touses)) > expelas = J(1, cols(shares), .) > } > if (quadratics == "") { > bofp = exp(lnp*beta') > lnpindex = lnp*alpha' > for(i=1; i<=rows(lnpindex); ++i) { > lnpindex[i] = lnpindex[i] + lnp[i,.]*gamma*lnp[i,.]' > } > } > else { > bofp = J(rows(lnp), 1, 1) // 1, so now div0 problem > lnpindex = J(rows(lnp), 1, 0) > } > if (ndemo > 0) { > cofp = J(rows(lnp), 1, 0) > for(i=1; i<=rows(lnp); ++i) { > cofp[i] = lnp[i,.]*(eta'*demo[i,.]') > } > cofp = exp(cofp) > mbar = 1 :+ demo*rho' > for(i=1; i<=rows(expelas); ++i) { > expelas[i,.] = 1 :+ 1:/shares[i,.]:* > (beta + demo[i,.]*eta + 2*lambda/bofp[i]/cofp[i]* > (lnexp[i]-ln(mbar[i])-lnpindex[i])) > } > } > else { > for(i=1; i<=rows(expelas); ++i) { > expelas[i,.] = 1 :+ 1:/shares[i,.]:* > (beta + 2*lambda/bofp[i]*(lnexp[i]-lnpindex[i])) > } > } > if (atmeans != "") { > st_matrix(expelass, expelas) > * > * For searches and help try: > * http://www.stata.com/help.cgi?search > * http://www.stata.com/support/faqs/resources/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/faqs/resources/statalist-faq/ * http://www.ats.ucla.edu/stat/stata/