program define ge_entr

    /************************************************************/
    /* computes some generalized entropy measures and bootstrap */
    /* confidence intervals                                     */
    /*                                                          */
    /* syntax: ge_entr <weightvariable> <incomevariable>        */
    /*                                                          */
    /* if no weights are needed, generate variable with all     */
    /* values set to 1 and use this as weight variable          */
    /*                                                          */
    /* for details see Biewen (2002): Bootstrap inference for   */
    /* inequality, mobility and poverty measurement, Journal    */
    /* of Econometrics, Vol. 108, pp. 317 – 342                 */
    /*                                                          */
    /* copyright Martin Biewen 2001                             */
    /************************************************************/

    version 6



    /* parse syntax */

    syntax varlist(min=2 max=2)



    /* set number of bootstrap replications and confidence level */

    tempname reps conf nobs

    scalar `reps' = 100

    scalar `conf' = 0.95

    scalar `nobs' = _N



    /* check if incomes > 0 */

    quietly sum `2'

    if r(min) <= 0 {

        di in red _newline "error: incomes must be strictly positive"

        exit
        }


    /* compute estimates */

    ge_ind `1' `2'


    tempname egem1 emld etheil ege2

    scalar `egem1' = r(gem1)

    scalar `emld' = r(mld)

    scalar `etheil' = r(theil)

    scalar `ege2' = r(ge2)



    /* preserve and initialize post */

    preserve

    tempname postnam

    tempfile boot

    quietly postfile `postnam' bgem1 bmld btheil bge2 using `boot', replace



    /* bootstrap loop */

    local i = 1

    while `i' <= `reps' {

        restore, preserve

        bsample

        ge_ind `1' `2'

        post `postnam' sqrt(`nobs') * (r(gem1) - `egem1') /*
                        */ sqrt(`nobs') * (r(mld) - `emld') /*
                        */ sqrt(`nobs') * (r(theil) - `etheil') /*
                        */ sqrt(`nobs') * (r(ge2) - `ege2')

        local i = `i' + 1

        }

    postclose `postnam'



    /* generate lower and upper confidence bounds */

    use `boot', clear


    sort bgem1

    tempname lbgem1 ubgem1

    scalar `lbgem1' = `egem1' - bgem1[`reps' - int((-(`conf'/2) + 0.5) * `reps')] / sqrt(`nobs')

    scalar `ubgem1' = `egem1' - bgem1[int((-(`conf'/2) + 0.5) * `reps')] / sqrt(`nobs')


    sort bmld

    tempname lbmld ubmld

    scalar `lbmld' = `emld' - bmld[`reps' - int((-(`conf'/2) + 0.5) * `reps')] / sqrt(`nobs')

    scalar `ubmld' = `emld' - bmld[int((-(`conf'/2) + 0.5) * `reps')] / sqrt(`nobs')


    sort btheil

    tempname lbtheil ubtheil

    scalar `lbtheil' = `etheil' - btheil[`reps' - int((-(`conf'/2) + 0.5) * `reps')] / sqrt(`nobs')

    scalar `ubtheil' = `etheil' - btheil[int((-(`conf'/2) + 0.5) * `reps')] / sqrt(`nobs')


    sort bge2

    tempname lbge2 ubge2

    scalar `lbge2' = `ege2' - bge2[`reps' - int((-(`conf'/2) + 0.5) * `reps')] / sqrt(`nobs')

    scalar `ubge2' = `ege2' - bge2[int((-(`conf'/2) + 0.5) * `reps')] / sqrt(`nobs')


    /* print out results */

    di _newline "----------------------------------------------------------------"
    di          "Generalized entropy measures with bootstrap confidence intervals"
    di          "----------------------------------------------------------------"
    di          "Confidence level = " %3.2f `conf' ",             " %4.0f `reps' " bootstrap replications
    di          "----------------------------------------------------------------"
    di          "GE^-1 = " in green %7.6f `egem1' in yellow  ",      confidence interval " in green "[" %7.6f `lbgem1' " ; " %7.6f `ubgem1' "]"
    di          "MLD   = " in green %7.6f `emld'  in yellow  ",      confidence interval " in green "[" %7.6f `lbmld' " ; " %7.6f `ubmld' "]"
    di          "THEIL = " in green %7.6f `etheil' in yellow ",      confidence interval " in green "[" %7.6f `lbtheil' " ; " %7.6f `ubtheil' "]"
    di          "GE^2  = " in green %7.6f `ege2'   in yellow ",      confidence interval " in green "[" %7.6f `lbge2' " ; " %7.6f `ubge2' "]"
    di          "----------------------------------------------------------------"





    restore

end



program define ge_ind, rclass


    /* returns r(gem1) r(mld) r(theil) r(ge2)    */


quietly {


    /* calculate products */

    tempvar w x wx w1dx wx2 wlogx wxlogx

    gen `w' = `1'

    gen `x' = `2'

    gen `wx' = `w' * `x'

    gen `w1dx' = `w' / `x'

    gen `wx2' = `w' * `x' * `x'

    gen `wlogx' = `w' * log(`x')

    gen `wxlogx' = `w' * `x' * log(`x')



    /* calculate moments */

    tempname m0 m1 mm1 m2 t01 t11

    sum `w'

    scalar `m0' = r(mean)

    sum `wx'

    scalar `m1' = r(mean)

    sum `w1dx'

    scalar `mm1' = r(mean)

    sum `wx2'

    scalar `m2' = r(mean)

    sum `wlogx'

    scalar `t01' = r(mean)

    sum `wxlogx'

    scalar `t11' = r(mean)



    /* calculate indices */

    tempname gem1 mld theil ge2

    scalar `gem1' = 0.5 * (`m1' * `mm1' /(`m0' * `m0') - 1)

    scalar `mld' = - `t01' /`m0' + log(`m1'/`m0')

    scalar `theil' = `t11' /`m1' - log(`m1'/`m0')

    scalar `ge2' = 0.5 * (`m0' * `m2' /(`m1' * `m1') -1)



    /* return results */

    return scalar gem1 = `gem1'

    return scalar mld = `mld'

    return scalar theil = `theil'

    return scalar ge2 = `ge2'

    }

end